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ABSTRACT 


The daveiopment of valid creep fracture criteria la a major topic In cur- 
rent engineering research. Two path-independent Integral parameters which show 
some diigree of promise as fracture criteria are the C* and (AT)^, integrals. 

The first portion of the present work reviews the mathematical aspects of 
these parameters. This Is accomplished by deriving generalized vector forms 
of the parameteri^ using conservation laws which are valid for arbitrary •three- 
dimensional, cracked bodies with crack surface tractions (or applied displace- 
ments) , body forces, inertial effects and large deformations. Two principal 
cunclusions are that 1-8 a valid crack-tip parameter during nonsteady as 

well as steady-state creep and that has an energy rate interpretation 

whereas C* does not. 

The dovolopment and appllcatAun of fracture criteria often involves the 
solution of boundary/ Initial value problems associated with deformation and 
stresses in etthci' laboratory specimens or actual components. Due to the power 
of the finite element method In treating complex geometries and non-linearities, 
it has often been used for this purpose. In the present work, an efficient, 
small displacement. Inf Initeslrocl strain, displacement based finite element 
model Is developed for general elastlc/plastlc material behavior. For the 
present numerical studies, this model Is specialized to two dimensional plane 
stress and plane strain and te power law creep constitutive relations. 

A mesh shift ing/remoshlng procedure is used for simulating crack growth; 

The model is Implemented with the quarter-point node technique and also with 
specially developed, conforming, crack-tip singularity elements which provide 
for ther strain singularity associated with the HRR crack-tip field. 
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GomparAaons are made with a variety of analytical solutlc-ns and alternate 
numerical solutions tor a number of problems. 
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INTRODUCTTON 

Fracture Prohloma and Fracture Criteria; A Review 
Characterisslng the dleplacemcnt, stress and strain fields associated 
with stationary and propagating cracks in solids characterised by various 
ideal ivLcd constitutive relations is one of the most important areas of study 
in fracture mechanics. The importance of those studies in not only that we can 
predict the stress or displacement fields in a cracked body, but also that 
knowing the nature of such fields we can possibly correlate observed fracture 
behavior with some aspect of these fields and thus arrive at valid fracture 
criteria. 

Cfitorift for Crack G rr , . Initiation 

The two macroscopic aspects of fracture for which correlations are coimr.only 
sought arc the initiation of crack growth and the rate of crack growth. The 
most notable liiltlatlon correlations are with the elastic stress intensity 
factor, K^, for the elastic (and/or small scale yielding) case [1,2] and with 
the J. -integral for cases in which plasticity may not be limited to the crack- 
tip [3,4]. The conditions under which those correlations are independent of* 

I geometry arc discussed in the cited references. The critical values of K, and 

i 

I for a given material ate denoted and respectively. It is implied 

by the use of the subscripts and "1" that these criteria are for the crack 
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opening mode (t.e., mode I). While similar criteria might be expected for the 
sliding and tearing modes (i.e., modes II and III), the experimental data for 
such studies is lacking. 

Both Kj, and have been shown to be crack-tip field parameters and both 
can be evaluated experimentally through energy considerations. The linear elas- 
tic, mode I crack-tip field determined in {5j shows that the asymptotic crack- 
tip fields arc proportional to K^.. Similarly, the asymptotic, mode I, crack- 
tip fields for power-law deformation theory plasticity have been shown (6.7] to 
depend upon the single parameters In the case of small scale yielding, Kj. 

is easily related to the energy release rate, (8], which is a measure of the 
potential energy decrease due to an increase in crack length. The quantity 
l>aa a similar potential energy interpretation in the case of deformation theory 
plasticity, and becomes identical to Gj. for small scale yielding of a stationary 
crack. 

In the foregoing discussion, the time dependence of the material's response 
and of the applied loading is assumed to be negligible. For creep crack growth 
these assumptions are no longer valid. We now consider crack growth initiation 
in materials which exhibit creep behavior, Whilo a significant number of creep 
fracture experiments have been reported in the literature, it appears that the 
primary Interest has been to find a creep crack growth rate criterion as opposed 
to an initiation criterion. As a result of this emphasis, many investigators 
use notched specimens rather than precracked specimens and many do not report 
data which could be useful in addressing the question of initiation. At 
present there seems , some indication (9,10} that when precracked specimens 
are used, the time required for creep crack growth initiation is negligible when 

^ The deformation theory of plasticity precludes elastic unloading from an elastic- 
plastic state and thus is mathematically equivalent to nonlinear elasticity. The 
crack-tip fields associated with power-law deformation plasticity are commonly 
referred to as HRR fields after the authors of references [6] and [7]. 
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cumparud to the il£e of the Gtpoclmen< It should be understood that this is not 
a universally acknowledged conclusion [11] and that further study is indicated. 

As noted previously, the second macroscopic aspect of fracture for which 
correlations with crack-tip fields are commonly sought is crack growth or propa- 
gation. The following suiranarlzes the development of criteria relating to this 
aspect . 

Criteria for Crack Propagation 

Slow crack growth occurring under constant load implies that the material 
response, is time dependent and is generally classified as creep crack growth. 

If the material's time dependent nature is negligible under the subject con- 
ditions, then it is assumed that crack growth requires an Increase in applied 
load. This latter case is typical of situations in which small scale yield 
conditions are not met and for which has been found to correlate with crack 
growth initiation. The primary Interest in this quasi-static mode of crack growth 
is that for some materials and geometries, the increase in load carrying capacity 
of the structure during quasl-statlc crack growth is significant. This implies 
that design procedures can be developed to take advantage of this added margin 
of safety. To justify such a procedure, however, there muast be some dependable 
means of predicting the crack growth versus load beliavior as well as predicting 
at what load the crack becomes unstable (t.e., is no longer quasi-static). As 
noted, is generally accepted as a valid Initiation criterion for this problem. 
For prediction of the subsequent growth, however, there are at least two proposed 
criteria which appear to provide reasonable correlations with experimental data. 
The first growth criterion can be stated as , where 

■ * 

J^^(Aa) is assumed to be a material property which depends on the amount of 
crack growth, a [12,13]. The subscript "R" denotes that this quantity 
characterizes the material's resistance to cracking. While strong theoretical 
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arguments can be given as to why this criterion should not be valid (except 
possibly for very limited amounts of crack growth 114)), it has been demon- 
strated that reasonable predictions can result from the use of this criterion 
for at least some classes of problems [12,15]. 

based on the theoretical objections to the use of except for limited 
crack growth, a second criterion which is based on the crack-tip deformation 
has been proposed [16,17]. This criterion results from finite element simu- 
lations of quasi-static crack growth experiments which Indicate that the crack- 
tip opening angle, CTOA (defined by the first finite element behind the crack- 
tip) becomes constant during crack growth. Whereas the CTOA, so defined, is 
clearly a mesh dependent quantity, the concept of ersek-tip defomatlons becoming 
constant with crack growth is physically meaningful. The procedure for applying 
this criterion in finite element based predictions of crack growth behavior is 
to use for initiation and for crack growth prediction until the computed 
CTOA has become constant with crack growth. Continued growth is then governed 
by this constant value of CTOA. Alternatively, a predetermined CTOA resistance 
curve can be used through'oUt growth. Crack growth Instability is assumed to oc- 
cur (for cither CTOA or as the criterion) when further increase in crack 
length results in the criterion for growth being exceeded without further increase 
in applied loading. The and CTOA criteria appear to provide reasonable 
correlation of ductile slow crack growth behavior for a variety of materials, 
geometries and load conditions [15,18,19]. 

Creep crack growth generally becomes a concern when components are operated 
at elevated temperatures. Whereas quasi-static crack growth can be on the order 
of mm/sec, typical creep crack growth rates are on the order of pm/sec. Compared 
to elastic-plastic quasi-static crack growth, the problem of creep crack initiation 
and growth is a relatively new area of study. 
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Numerous exporlmcntal studies have been undertaken with the purpose of 
findirg a parameter which correlates with creep crack propagation rate* (See, 
tor example, the review article [20] and [21-24].) Most of these investigations 
consider as candidate parameters, K^, some form of net section (or reference) 
stress, and In more recent studies C^. The C* parameter Is the steady-state 
creep analogue of (In the sense of [25]) in that the definition of C* is the 

same as that for except displacements and strains are replaced by their res- 
pective rates [2Cs). 

It is illustrated In Fig. 1.1 that the above three parameters can be expected 
to correlate three distinctly different creep crack growth situations. In Fig. 
1.1a, a crack and its associated ligament are shown for a material and geometry 
which results in negligible creep strains everywhere except in the vicinity of 
the crack-tip. This condition is analogous to that of small scale yielding in 
elastic-plastic fracture. Fig. 1.1b represents a situation in which C* might 
be considered an appropriate parameter. This situation is characterized (i) by 
the body being essentially at stesdy-stace .’’•eep conditions (which implies very 
slow crack propagation) and (ii) by the creep-damage process-zone being local 
to, and therefore controlled by, the crack-tip field. Pig. l.lc Illustrates 
the type of situation for which net section stress might be expected to control 
crack growth. In this case, the main feature Is the widespread creep damage 
zone. 

It is seen from Fig. 1.1 that Intermediate situations can occur. For 
example, suppose a particular material and geometry results in a crack -'ropagation 
rate such that elastic strain rates are not negligible compared to creep strain 
rates (l<e., nonsteady creep) and at the same time, creep strains are no longer 
localized to the crack-tip region, l^ile neither or C* could be valid 
parameters for this case, it appears reasonable to expect that crack growth 
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0. Kj controlled 
behavior 



b. C ,(aT^ controlled 
behovlor 


c. a net controlled 
behavior 


region with 
creep damage 


Fig. 1.1 


Conditions for which creep crack growth 
parameters are expected to be valid 


r»Ce Is still deterrolnod by the local crack-tip field since the creep domsge 
process zone is still assumed to be local to the crack-tip. 

A paramatcr which is apparently capable of spanning the gap between 
controlled grtv'xh and controlled growth has been introduced [27]. This 
parameter is referred to as (AT)^ and is defined by a path-independent vector 
Integral. A detailed discussion of a generalized (i.e., C*) and (AT)^ 
is given in Section II of this work with a principal result being that the 
energy relationship commonly used for experimental measurement of C* does 
not apply to C* but rather applies to the (AT)^ parameter. This means that 
the experimental results are actually showing a correlation with (AT)^ rather 
than with C’*'. Based on the theoretical validity of (AT)^ as a crack-tip field 
parameter for nonsteady as well as steady-stote creep and based on the mounting 
experimental evidence thot crack propagation rate correlates well with (AT)^, 
it seems the creep crack growth rate problem la close to having a solution. 

Motivation for the Present Work 

In the following, we review previous studies to the extent required to 
place the present study in perspective and briefly introduce the present work. 

The nonlinear nature of creep constitutive relations precludes analytical solu- 
tions for either stationary or propagating cracks in a creeping material. For 
stationary cracks in a power-law creep material, however, it is known that 
the HRR fields are present in the vicinity of the crack-tip [26]. (Since the 
singularity in creep strain rates is greater than that in the stresses, and thus 
elastic strains, it follows that the HRR field exists at the crack-tip during 
nonsteady as well as steady-state creep.) For propagating cracks , it appears that 
the HRR fields no longer exist at the crack-tip, but that analytical tools exist to determine 
the fields which do exist [28,]9|. Wliilc knowledge of the crack-tip field is 
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valuable, the solution of boundary value problems mast depend on numerlcsX 
methods. The finite element method, in particular, shows promise for solving 
creep crack growth problems. 

Only a few studies on finite element modeling of creep crack growth have 
been reported. The earliest is apparently that of Chtani and Nakamura [30]. 
This study simulated crack growth with a node-release technique and assumed a 
criticol crack-tip plastic strain criterion for creep crack growth. The rate 
constitutive l»,v contained an elastic term and a creep term based on the gen- 
eralization of the cnlaxlal Norton power law, 

Hlnnerchs [31] uses the Bodner-Partom constitutive law [32] and a node- 
release technique for modeling crack growth. In this work, sev^^ral candidate 
criteria are examined by simulating crack growth experiments. Due to the ap- 
parently limited crack growth (<0.5 mm), the short test durations (one hour) 
and the lack of crack growth measurement data (which requires the development 
of a so-called h'/hrid experimental-numerical procedure to estimate the crack 
growth history) , it seems the general applicability of the conclusions from 
this study are questionable. It also seems likely that the methods for evalu- 
ating C’'*' iu this study are incorrect^ and thus the conclusions concerning C* 
should be regarded accordingly. 

Ehlers and Riedel (33] have conducted a finite element analysis of a 
stationary crack in a compact specimen. The primary emphasis in this study 
is on the nature of the crack-tip field during the transition from the Initial 
elastic field to the steady-state creep field. 




Wliile the details of the numerical procedures for evaluating C* are not given 
In (31), it appears that the W'^ term of C* (see Section II) is incorrectly in- 
terpreted as a history dependent quantity as opposed to a quantity dependent 
solely on the steady-state stress and strain rate. 
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The strength of the HRR field during the transition period Is determined 

through fitting the iicnr tip equivalent stress field. The calculations use 

eight~aodcd isoparametric elements with quarter-point elements being used at the 

- 1/2 

crack-tip so as to have an r strain singularity. Creep crack growth and 
creep crack growth criteria are not considered in this study. 

Thu finite element equations for the creep crack growth model being used in the 
present study are derived from the principle of' virtual work in Section III. 

Section V presents the results of several analytes Involving both stationary cracks 
and propagating cracks. The creep crack growth simulation Is via a mesh shifting/ 
reraeshing procedure. Calculations are made using the quarter-point element tech- 
ni<]ue as well as with a specially developed (Section TV) compatible clement which 
incorporates the HRR, j strain singularity. 

An important aspect of the current work is the study of the (AT)^ parameter. 

In particular, the meaning of (AT) , its relationship to G*, and its calcula- 
tion within the context of finite element analysis nr explored in depth. 

A series of crack propagation calculations are combined with analytical and 
experimental results in Section V to show that creep crack growth in 304 stain- 
less steel at 650°C occurs under essentially steady-state creep conditions. This 
implies that the crack growth rate for a given crack length and load can be 
determined from a steady-state creep solution which does not depend on the 
previous load and crack growth histories. This observation implies that 
simple crack growth prediction methodologies may be developed. 
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SECTION II 


DERIVATION OF THE (AT) , J AND C* INTEGRALS 

****** C ***** 


Preliminaries 

We shall consider problems which exhibit the following constitutive be- 
ll av lor : 




( 2 . 1 ) 


We denote the cartesian coordinates of the undeformed body as x^. Defining u^ 
as the rate of displacement (or velocity) of a material particle from the cur- 
rent configuration, then e , is the symmetric part of the rate of displacement 

T ’“l 

gradient s (^“) ^ "*■ • The gradient operator is with res- 

pucu to the current coordinates y^, where it is understood that 

Lljki tensor of instantaneous elastic moduli. We let denote the 


corotatlonal rate (or ’’Zareraba-Jaumann rate") of the Klrchhoff stress CJ 


9y, 


ij 


where o.^ is related to the Cauchy stress t . by o ■ JT, . (j»det [•--”]) . The 


U 


IJ 


Ij 


3x 


equivalent Klrchhoff stress o is related to the devlatoric Klrchhoff stress 

©q 


"Ij '■ "iJ - 

those of the familiar Norton's law 


1/2 


The parameters Y and n are 


e • y(o ) 
eq eq 


were 


s ,- 


1/2 
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Wo will uae the notation; denotes a second order tensor; (_) iBiplies 

a vector; ^ B • £ implies a >• B ,c ; A •• B . G implies A., • B..C. . ; 

A:B “ A,,B,,. Also note that ^B;C Implies -r-^ C and V . B implies 
ij ij -t- - ' 9y^ jk -t 


5B 
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A Conservation Law for Finite Elastic and 
Nonsteady Creep Material Behavior 

The discovery of conservation laws and the possibility of deriving path- 
Inciupendent integrals from these laws are not particularly recent occurrences as 
discussed in [34]. However, the literature in this area has been rather piece- 
meal and tlierofore difficult to assimilate. The recent work of Atlurl [27) has 
done much to unify and generalize this subject and is the basis for the following 
presentation. 

We will consider a very general conservation law which has been given by 
Atluri, but will limit our discussion of this law to materials characterized 
by (2.1). We will use cartesian coordinates exclusively. Note that by 
special selection of material constants (i.e., Y’*0) , (2.1) can be specialized 
to elasticity. Alternatively, by assuming that the stresses are invariant with 
time, (2.1) can be specialized to steady-state creep behavior. 

In the following presentation, the current configurat-ton (i.e., the config- 
uration at time t) is the reference configuration. There may be initial stresses 
existing for this reference configuration. If stresses do exist, then they are 
assumed to satisfy the linear and angular momentum balance condition (i.e., 
equilibrium) 


• I + " 05 I " 

where p^,, ^ and ^ are the current mass density, body force vector and accel- 
eration vector. 
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A consecvatlon integral relation given by Atluri [27] for a closed volume 
(at the cur, 'rent time, t) , which is free from singularities and any other 
deCects (which would preclude the application of the divergence theorem), is: 

0 » f {V AW - (^T):Ag ^ V . f(x+At) . Ae) (2.2) 

J\j t ^ ^ 

t 

- Pj.(f.-a) • ^e)dV r [n^ . (j+At) - £] . AadS 

Js^ 

f Su. • (T+^t) • (Ae-Ae)dS 

G 

In (2.2), At is the Incremental f Irst-Piola-Kirchhoff (nonsymmetrlc) stress 
(At ■ [Ag - Ae . g]/J), where Ag is the material increment of Klrchhoff stress. 
The current mass density is denoted p^, and jE and a are the body force and 
acceleration vectors at time t+At, respectively. and are the portions of 
the boundary of upon which prescribed tractions, E, are acting and at which 

pi-escribed displacement gradients, Ae, exist, respectively. The current outward 
normal to S or S is n . The quantity AW, discussed in detail in [27], is the 
Incremental stress-working density in time At, and is given by: 

1 T 

AW “ x:Ae + jAt :Ae = x:Ae + AU (2.3) 

where 

AU ■ •|•At'*':Ae (2.4) 

The validity of (2.2) is readily verified through the two identities [27]; 

V. AW - V (x:Ae) + _^A0 = ^^x;e (2.5) 

4- ^Ae;x ^Ae:At’^ 

“t ' - — t - - 

and 

,ij. • [<I+At) . Ae] - [_^ . (x+At)] . Ae 4- T^Ae: (x+At)^ (2.6) 
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the satisfaction of linear momentum balance In V : 

t 

ij. • “ 0 (2.7) 

and the satisfaction of the boundary conditions:^ 

Jit • " i on (2.8) 

Ae " Ae on (2.9) 

Note that Identity (2.5) assumes that T ( the initial stress for the Incre- 
ment) Is an explicit function of its position in V^, The existence of AU 
is shovm and discussed in the work of Atlurl [35]. 

Having the relation (2.2) It Is now possible to specialize this relation 
to finite elastic behavior or to steady-state creep behavior. However, since 
we are primarily interested in the path- independent integrals which can be ob- 
tained from (2.2) we will postpone the specialization till after we have de- 
rived the general path- independent integral (AT) . 

— c 

Path-Independent Integrals for Fracture Analysis 

The conservation integral (2.2) is used [27] to obtain a path-independent 

integral which is applicable to the analysis of cracks by considering a volume 

V - V such as illustrated in Fig. 2.1. (Note that a two-dimensional case 
t 6 

is illustrated for simplicity). The use of the divergence theorem for the 
region depicted in Fig. 2.1 results in (2.2) being rewritten 



. (x+At) 


-P (f-a) 
t 


Ae]dS 


Ae]dV 


( 2 . 10 ) 


The validity of (2.2) does not require where dV denotes the surface 

bounding Vj.. Therefore, dV need not coincide with the boundary of the body 
under consideration. 


13 - 


/ n AWdS + / n. AWdS - / 

j. -t Ji> Js 

12 45 ( 

I EL • (f+At) . AedS 
js -t 


AedS 


“ I fa 

Jr 


_ AW - n . (x+At) . AeldS 2 (AT) 
— t - - ~ — c 


( 2 . 10 ) 


In welting (2.10) it has been assumed that implies that 


‘234 


docs not coincide with any exterior boundaries. This has been assumed 


purely for convenience of notation. We have also used the notation 
Noting that (2.10) contains two equalities, it can be verified by inspection 
that (AT)^ dep^ds on £ (or more generally F^) but that it does not depend 
on the selection of F^^^* In this sense (AT)^ is path- independent (i.e., 
Independent of the selected £ar-fleld path). Following the reasoning of Atluri 
[27], we define (AT) as the limit of (AT)*^ as e goes to zero.^ 




I 


AW - n . (l+At) . Ae]dS 


( 2 . 11 ) 




[l> AW - n , (t+Ac) . AeldS 


234 

Lt 
€“>0 


I f [(-V x):Ae - P,.(jL-a.) • Ae]dV 

( Jy -V -t- - t 

' t € 


f n AWdS + f n AWdS - / E . 


AedS 




(x+At) . AedS 


I 


The existence of the limit is shown in Appendix A. 
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Fig. 2 




1 Contours for applying the conservation law to 
a two-dimensional, cracked body 
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Dy defining (AT) as the limit of (AT) as € goes to zero, it is seen that 
(AT)^ is entirely determined by the asymptotic near-tip fields. It will be 
shown later that the converse to this statement is also valid when the near- 
tip fields are the HKR fields. That is, it will be shown that (AT) entirely 

•MM ^ 

determines the asymptotic near-tip fields. 

Often it happens that only the first component of the vector quantity 
is of interest. We will write the first component of (AT) as (AT.) . Also, 
the quantities (T) and (T, ) will often be used in place of (AT) and (AT.) t. 
These quantities are related by 


(T) 


c 


Lt 

At ->-0 


(AT)^ 
~ c 

At 


( 2 . 12 ) 


However, in the presentation of numerical evaluations of (AT.) we use (T, ) 

X C X c 

as a convenient, approximate notation for (ATj^)^Mt. 

We now consider two special cases of (2.11). For symmetrical deformation 
about the axis and cracks oriented along the x^ axis with traction free 
crack surfaces, no body forces and negligible inertial effects, the first 
component of (AX)^ is 




Lt 

">0 



[nj^AW - nj(T^^+Atj^)Aej^j^]dS 


(2.13) 



(nj^AW - nj (tj ^+Atj j^) Ae^^j^ldS 


3t 


iiAe. 

1 


9yi “'ji 


dV 


Note that the limit of the volume integral has been written in its explicit* 
form as a result of the arguments for the existence of this limit, given in 
Appendix A. If, in addition to the above conditions, the strains are in- 
finitesimal and the deformations small, then there is no need to differentiate 


betwomi and y^, becomes identical to we have (recalling 

Ac^j S 3Au^^/9y^)s 


It r 

''i>c-.-o 1V“ - 


(AT 


(2.14) 


X 

' 234 

-X. 


8Au 

IniAW - nj(T^^4AT^^) 


OdS 


ati. 

-r — Ae , dV 

3x, ij 


The replacement of volume integral is made possible by 

the symmetry of and thus ax^j/^xj. 

Physical Interpretation of (A^)^ 

It has been shown by Atlurl [27], that the vector (AT)^ has the following 
physical meaning. Let two bodies with non-propagating^ cracks be identical 
except for the second body having an additional, arbitrarily directed, in- 
finitesimal increment in crack length characterized by the vector It is 

assumed chat both bodies experience identical load histories. Define total 

2 

potential energy increments corresponding to the time Increment At as 


AEj^ = Ai|i^ + A«j^ -f AKj^ 


AE^ = A»|<2 + An^ + AK^ 


(2.15a) 

(2.15b) 


for the first and second bodies, respectively. In (2.15), -Aij) is the incre- 
mental work of external forces. Ail is the Incremental stress-work and AK is 
the increment in the kinetic energy. (It slx^uld be noted that AJl includes the 
inelnstically dissipated energy.) Then 


^Atlurl [27] has shovm that the l/r singularity in kinetic energy, which is as- 
sociated with dynamically propagating cracks, changes the interpretation. 

2 

Note that sign convention for AEj^ and AE^ is opposite to [27] so as to conform 
to conventional usage. 
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(AT ) bdc, 
1 c J 


-(AE^ - AE^) 


(2.16) 


where b is the length of the crack front. 

If one Is only interested in self-similar crack extension in the x^^-directlon, 


then dc. 


dc. 




0 and 


AE^ - AEj^ 
bdc. 


(2.17) 


Therefore, (AT)^ is relai.:!d to the incremental potential energy difference 
between two bodies which are identical except for an Incremental crack length 
difference dc . 

Finite Elasticity and ^ 


As noted previously, the constitutive law (2.1) can be specialized to 
elastic behavior by choosing y to be zero. Therefore, (AX)j, o® defined by 
(2.11) is a valid crack-tip characterizing parameter for general nonlinear 
elasticity with finite strains, large deformation, body forces and inertial 
effects. Howerver, the basic premise of clastic behavior is that the con- 
stitutive relations are independent of the histories of deformation and stress. 
TIUs means that the constitutive relations can be derived from a potential. 

For instance, a potential, U, exists for jt, the first Piola-Kirchhof f stress, 
such that 


T 8U 

t m 

5 3e 


(2.18) 


In the following, we consider the reference configuration to be the stress-free, 
undeformed configuration at t**0, and therefore drop the subscript t for con- 
venience. As a result of the existence of the relation (2.18), it is possible to 
state two identities which are analogous to those of (2.5) and (2.6) for the rase 


tSO. 


VU - 


DU 


De 


mn 


De Dx, 
mn 1 


7e;t 
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(2.19) 


and 


V . (t . e] - V . t.e + Vc;t (2.20) 

Similarly, the linear momentum balance (l.e., equilibrium) condition is now 

V . r, + P(f. - a) - 0 (2.21) 

and the boundary conditions 

n . t ■ t on 


e ■ e on ^ 


( 2 . 22 ) 

(2.23) 


Noting the similarity o£ equations (2.18) through (2*23) with (2.4) through 
(2.9) it IS easy to arrive at the following conservation law 


0 - / {VU - V . (t , e) - p(f - a) . e)dV 


(2.24) 


L -- 


. t - It ] . edS + 


Is"- 


. t . (e “ e)dS 


Following the procedure used 'in deriving (A^) from (2.2) we apply the diver- 
gencc theorem to (2.24) for the volume V - and take the linlt es t goes 
to zero. The path-independent vector quantity resulting from this procedure will 
be called J. 


J 



in . t . e]dS 


(2.25) 



In writing (2.22) through (2.25) it is understood that o is the mass density 
in the reference configuration, n is the unit normal in the reference configu- 
ration, £ are arbitrary body forces per unit mass, a is the absolute material 

- n - 


acceleration, and ^ arc prescribed tractions per unit undeformed area S . 

We now consider several special cases of (2.25). If the problem being 
considered Involves a crack oriented along the coordinate direction and Is 
loaded so that only mode I crack-tip behavior occurs, then la of primary 
interest and we have: 


•’i ■ X '"i" ■ "I'ij'ji''"'' * ‘“I Xv „ ‘“i - 

- X ■ X 

t e I 

If in addition the problem involves infinitesimal strains, small displacements 
and traction-free crack surfaces, we have 



(2.27) 


where use has been made of the existence arguments of Appenix A in taking the 


limit of the volume integral.. 

For elastic behavior and non-propagating cracks, Atlv.rl [27] shows that 
J has the meaning of energy release rate to a process zone in the sense that 


J, bdc, *= 
k k 



dt 


(2.28) 


where b is the length of the crack front, 

DE Dll) Dn DK 

!- a SL ^ _ — i. ^ L. 

Dt Dt Dt Dt 


and 

I dt - E (t + dt) - E (t) 

c c 

For an elasto-dynamlcally propagating crack (i.e., singular kinetic energy) 
Atluri [27] concludes that 



(2.29) 


J, bdc, 
k k 


■[iU 


wliere L is the Lagrangian (l.e. , L 5 K - - Q ) such that 

( c c c 


DL. 

f 

"dT 


D<p 

t 

Dt 


DK^ 

"dT "Dr 


TherofOLC, J|^ has the meaning of "rate of change of Lagrangian per unit crack 
growth". 

We now consider the special case of Bteady>state creep behavior. 

Steady-State Creep and C* 

It hae been shown that (AT) characterizes the crack-tlp field for mat- 

— c 

orlals which exhibit creep behavior such as In (2.1). it is known that under 
cectain conditions of applied loading, the constitutive relation (2.1) can (after 
long times) result in a steady-state. This steady-state is primarily charctcr- 
laed by the time independence of the stresses (l.e., AU • ■ 0), Specializing 

(2.13) to steady-state conditions, we define the steady-state value of (^I.j^)j,5 

(2.30) 






■L 




- r 
L 


dV 


*234 t 

Because (2.1) results in a power-law relation at steady-state, which is 
analogoufj to the power law deformat Ion- tfpory plasticity (or essentially non- 
linear elasticity), Goldman and Hutchinson (26) have suggested a path- Independent 


C| Integral, 


r r 

1 Jp [ 1 "j^lj 3x^ 


dS 


(2.31) 


where 


W* 


T dc 

Jo 


(2.32) 
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The question of how C* and (AT ) are related, la a natural one. Before 

JL X CSS 

obtaining an equation relating C* to (AT,)_^ , however, the conservation integral 

X X css 

(2.J) will be used to derive a generalized vector integral jC*. 

In specializing (2.2) to steady-state we note that now stress is a 
function of the strain rate, and that stress increments arc zero. Thus, 

AW t:Ac. Also we may write: 


£ 


[7 AW - (rT):Ae)dV - / t: 

\ ^ -- - L - 


7 AedV 


Thus, at steady-state, we may write (2.2) as: 


0 ■ I tls(2*Ae) - ^ . (r.Ae) - p (f - a) . A§)dV 
J -t -t - « t - - 


/ (n . T - t] . AedS + I n.. . T . (Ae - Ae)dS 

S ~ 

t e 


e)dS 


or equivalently, in rate form, 

0 - y* tT:(V^e) - ^ . (T.e) - P^(f - a) . ejdV 
^t 

f tlL. • X - t] . edS + / n . T . (e - e 

t e 

Using the symmetry of T v note that: 
x:7^e » x:^r (a + e^)l ■ x:7 c 

~ — t ' ~ — t ^ 2 “ - ml •« -i: 


As a result of the incompressibility condition (e^^ '■ 0) we have 
t: 7 e ■ T ' : V c 


and 


i 


X' d^ 


Ij "ij 


(2.33) 


(2.34) 


(2.35) 

(2.36) 


which leads to 


(2.37) 



Combining the results of the above manipulations we have: 


„ ... _ m m " mn „ . _ „ . 

V w* = -r — • -tT r - I “ T:7^e. 

-t 3y. 3e 3y. *' -t- ~ — t- 

•'l mn 'I 


Using (2.38) and the divergence theorem while applying (2.33) to - 
define the vector quantity (£*)^; 

/ tn.W* - n. . T .e]dS - f P (f - a) . edV 
Jp — t . ~ Jy _V 

234 ''t € 

+ I n W*dS + j n W*dS - / it . edS - / n . t . IdS 

■^r,. Jr,. Js. ” ■ Js - 


Vj. we 


fr 


n Wi*' - n . T . e]ds s (C*) 
— t — t ~ — 


If we define the limit of (£*) as e-K) to be C^*, we have a quantity which 
characterizes the crack-tip field and is Independent of the selection of ^ 23 /i' 
Restricting our attention to problems involving symmetric deformations about 
the axis and cracks oriented along the axis, with traction-free crack 
faces, no body forces and negligible inertia effects, we find that 


C* = 

1 C-+-0 


I 




(2.40) 


. f 

J r 


- "j 


In computing W* it is convenient to invert (2.1), substitute the re- 
sult into (2.36) and use the following identity to complete the integration: 


14-n 

d(e ) n 


o 1 . 1-n 

T (' ) (€ ) n 

3 n eq ij 


(2.41) 



Tha result of this manipulation is: 

1 


Relationship of (T ) and C* for Steady-State Creep 

Jw C 8 8 A 


l+n 


(2.42) 

(2.43) 


Now we will relate CJ of (2.40) to the steady-state value of 


First wc rewrite (2.30) in rate form as: 


(X ) « 

^ I'css c -^0 


X ‘"I'u'u - ■’ 




(2„44) 


" f '*‘ 44^44 “ "»''^4 4^41 “ f ■ e dV 

•Jp ‘ 1 ij ij j ji il Jy 3y^ ij 


234 


Using the notation W“T,.e.. we have 

ij ij 


(T.) 

1 CSS 


■ X '"1“ ■ "j'li'iU 


ds 


234 
9t , 




dV 


(2.45) 


Noting that; 


“ " "ij'u ■ ’u 2 ‘'u ■ Uj'13 


(2.46) 

it is seen that W is the rate of stress-working density , while W* is just a 


mathematical potential for t| . As a result of incompressibility we can write 

J 

I 1+n 


W •• 0 c ■ y(o ) _ (— )” (e ) "^ 

eq eq eq Y eq 


(2.47) 

as contrasted to W* of (2.42) and (2.43). Comparing the left equalities of 


‘^Thls result is only valid for steady-state creep and is obtained through the 
Substicution of the steady-state specialization of (2.1) into (2.46). 
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(;i. 'tO) ami it Is seen that (T ) and C* are related by: 

1. C S3 S 4' 


• 1 ^ 
(T.) - C* 

1 OSS 1 t ^0 


f "i<“ - 


W)dS 


(2.48a) 


""ir 'fn f )^^^dS 

1 n+i t >0 Jp 1 eq 
e* 


(2.48b) 


Appomlix IJ gives several numerical examples of relation (2.48) for two rather 
extreme vaUiea of n. 

It 1.8 uow clear that jG'*' and (A^) are not equivalent quantities under 
any ooiulltion despite their being derivable from the same conservation law^ . The 
quantity (AX)^, follows more directly from the conservation law and is the more, gen- 
eral quantity not only in tha^; it t.s applicable to nonsteady us well as steady- 
utute creep but also in that it is applicable to constitutive laws which are more gen- 
eral tlian (2.1), The quantity relies on the special property of (2.1) which allows 
the exl.stencu of a potential W*' for the stresses (t*). Furthermore, since W* does 
not have any physical meaning, whereas W has the meaning of stress-working density, 
it la uiulorstiindablc that (AX)^, energy interpretation whereas C''' does not. 

It is for this reason that it seems more appropriate to refer to experimental measure- 
• 

cl ^ * 

ru'uts of - as measurements of (1\) as opposed to measurements of Gi' or J . 

da 1 c 1 1 

The URR Field 

We now give the HRR field in terras of (AT,) . Whereas similar relations 

1 c 

have been written in terms of C* for steady-state creep [36], the relations in 

terms of (AT^)^ will be valid for nonsteady creep as well as steady-state creep. 

The HRR field as given in [37] but modified for creep by replacing c and u by 

tj ^ 


and respectively, is: 


Note that these equations arc derived on the assumption that tij^O (i.e., creep 
steady-state). Therefore, in order to have a well defined creep constitutive 
law we must have Yf*0 ^f'd n finite. 
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(2.49o) 




n+1 


|ay(o).a^^(o)i 



YK r 

f. 


-n 

n-t-i. 



(2./*9b) 


u “ YK r 

X ^ 


n+^f 




(2.49c) 


where 3 (0) has been nbnnallzed to have a maximum value of unity and K, and 

eq •' o 

K are amplitude factors which are related by 




(2.49d> 


It cun tlierefore be seen that the asymptotic crack-tip fields are entirely 
determined when (or K^) is known or specified. Combining (2.49a,b,c) with 
the first equality of (2.14), using (2.49d) and rearranging, gives: 



(2.50) 


where is analogous to I defined by Eq. (24) of [6] except for the factor 
u/ (n+1) multiplying the energy density term. To be explicit, 

I* - I + -rr r (ft)]"’’’^ cosOdO (2,51) 

n+1 eq 


It is therefore seen that knowing the value of (fj^)^ is equivalent to knowing 
K and thus is sufficient for defining all aspects of the asymptotic crack-tip 
field during nonsteady creep as well as under steady-state creep conditions. 
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SECTION III 


DERIVATION OF FINITE ELEMENT EQUATIONS 


Before stating the variational principle and deriving the finite element 
equations, it is beneficial to Illustrate the constitutive behavior to which 
the finite element model is addressed. 

Elastic/Viscoplastlc Constitutive Relations 

A rather general rate constitutive law proposed by Perzyna [38], can be 
written in incremental form as: 


where M, E and v are the elastic shear modulus. Young's modulus and Poisson's 
ratio, respectively, ^ deviator Ic stress and y is a 

viscosity constant of the material. In writing (3.1) it is Implied that 


Ae 


ij 


e , 

Afiij + 


vp 

ACij 


where At\. and Ae^^ are the elastic and viscoplastic strain Increments, rea- 
ij Ij 

pectivaly. The yield function governs the magnitude of the incremental 

viscoplastic strains through the function <4i(F)> where 

0 for F ^ 0 
't(F) for F > 0 


«t(F)> 


(3.2) 
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The relative magnitude!} for the incremental viscoplastic, strain components 
are seen to depend on the factor This implies that when con- 

sidered as a vector in nine-dimensional stress space, is always directed along 

the normal to the surface iC'T , . . 

ij k 

By choosing F t f « “ (3/2(t such that 3f/3t^j - (3/2) 

(I'y/o^q), and choosing 4>(F) “ (F)” we find that (3.1) becomes 


Ac 


ij 1 A • ^ 1-2V 

2m “^ij ■ 3E 




(3.3) 


This represents the special case of creep behavior which is considered exclusively 
2 

in this study. It should be understeod, however, that the finite element 
model which is described below is applicable to the more general behavior rep- 
resented by (3.1). 

De rivation of Finite Element Equations 

The finite element model is derived from the principle of virtual work 


/ T,,6€,,dV - f t 6u dS 

Jv ij Js i i 


0 


(3. A) 


IM the present finite element analysis, we assume only Infinitesimal defor- 
mations and strains; hence there is no need to dif f erentiate between the deformed 
and undeformed configurations. 

In writing (3}.4) it should be noted that x.. are the stresses 

ij 

existing at time t 4* At (where t is the current time), t^ are the prescribed 

tractions on S at t 4* At , and 5n [ 6€ « (1/2) ((5u .4-6u. .)] are arbitrary 

0 i ij l,j j,i 

compatible virtual displacements. 

Following customary procedures we introduce the element displacement 
shape functions which relate element displacement u^ to element nodal 


Xt is common to choose f^F in which case we have what is called an associative 
law. 

2 

While (3.3) Is known to deviate from real material behavior (especially for primary 
creep) it is a widely used constitutive law and therefore has been adopted in the 
present study. 
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diaplacemonts Iq} 

u “ (u) ■■ tN]{q}; 5u “ {<Su) >■ lNl{*5q) (3.5) 

1 i 

We again use the customary notation wherein strain (and stress) components 
are placed In one-dlmenslonal arrays 

{<) - fBHq); {&€} - [B]{6q) (3.6) 

Substituting (3.5) and (3.6) into (3.4) and applying conventional procedures 
f.oK assembling element matrices Into global matrices we have 

l«L[X 

- f tu'^[N]ds]jj - 0 

s 

0 

e 

Since {iSq} are arbitrary virtual nodal displacements, it follows that 


ILL 

L {tl'^'lBldV - 

f {t}'^[N]ds]| 2 {jf - {0}^ 

(3.7) 


‘''^e 

e 



We now express the stresses {x} at t + At in terms of the current 
stresses, and the Incremental stresses corresponding to the time in- 

crement At; 


In (3.8) and in the following, the I and l+l subscripts designate the 
incremental solution with which the quantity is associated. Application of 


the incremental elastic constitutive law results in 


« {t}^ -f [E]({Ae}j^_^j^ - 


vp'1+1 


(3.9) 


where (At } , are the incremental viscoplastic strains and (E] is the matrix 

vp ItI 


of clastic constants. Substituting (3.9) Into (3.7), taking the transpose, 
and placing the known terms on the right hand side we have the final form of 
the finite element equations: 

^^vp^I+1 ■" (3.10) 

where 


IKI 


ele 


f [B]'^tE][BldV 
Jv 


(T) 


t 

1+1 gIu 


j [N]'^{t}j^^dS 


}txi ~ f f [Bl^lEHAe 

vp I+l ele vp I+l 


(3.11) 

(3.12) 


(3.13) 


{R}, 


s r 
Jv 


[B]MT}j,dV 


(3.14) 


The above volume integrals are evaluated in the current work by 2x2 Gauss 
quadrature. The array CT) is input directly in terms of node point forces. 
Solution Procedures 

It should be noted that [K] of (3.10) is just the elastic stiffness 
and therefore only needs to be formed and decomposed^ once. This results in 
significant savings in the number of computations per time step as compared to 
methods using stiffness matrices which must be reformed at each step (i.e., 
tangent stiffness methods). It sliould also be noted that the term {S },,, 

Vp XT'! 

is computed from incremental vlscoplastic strains which are esti - 

mated Using {t}_ and the material constitutive law (3.1). Only for the 


The equations (3.10) are solved in the current work by the decomposition [K]» 
[L][D][L1^ where [D] is a diagonal matrix (the only nonzero entries are those 
on its diagor.al) and [L] is a lower triangular matrix (the only nonzero entries 
are those bejow its diagonal); see for example [39]. 
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siH'Cial Mituation when the stresses do not change with time will this estimate 


bo exact. Having obtained the Increraental nodal displacements by 

eolvinK, (3. 10), one can easily find the total Incremental strains {Ac)„., via 

1+1 

the Incremental analogue of (3.6) . We now describe two procedures for ob"- 


tainlng Ir} 


1+r 


The first and simpler method to obtain is to substitute the estl- 

I+l 

mated used in solving for into (3.9). If one does this, then 


it happens that 


^^^ 1+1 ^*^’^ 1+1 


(3.15) 


(3.16) 


and therefore (3.10) becomes for the next step^ 
lK]{aQ)j.^2 " ^^^1+2 ^^vp^I+2 " ^^^1+1 

” ^^ p ^+2 

This method was compared to the following method and was found to require smaller 
time steps to achieve similar results. 

Rather than using the estimated values of {Ac and (3.19) to 

detormlne {t}j.^ 2> constitutive relation (2.1) is integrated over the cur- 
rent time step at each Gaussian quadrature point with the condition chat total 
strain It) varies linearly with respect Co time from (t) to (f) .. (The 

X. 4* » X 

present study uses an Eulerian scheme with each time step being divided into 
five subincrements.) The result of this procedure is better adherence to the 
pustulated constitutive law at the expense of introducing a somewhat unequili- 
brated stress state. The amount of disequilibrium depends on the accuracy of 
the original estimate for the Incremental viscoplastic strains and thus on the 
time step size. ’ 

At this point one has two alternatives. The first is to use the 


4 


This procedure results in Che current model reducing to that of Zienkiewicz 
and Corraeau [40], 



viscoplas^ Ic strain increments obtained through the time integration procedure 
as an improved estimate and to re-soive (3.10) for the current time step. This 
procedure would, after several iterations, result in a stress state which is 
equilibrated to within some small user specified tolerance. With this type of 
procedure the time steps could be as large as those used witli tangent stiffness 
methods. Further, it is reasonable to expect the solution to be at least as 
accurate as if a tangent stiffness method were used. 

The second alternative is to go Immediately to the next time step with 
the understanding that the term {R}j. in (3.10) results in the disequilibrium 
from the T.th step being corrected in the 1+1 step. This feature is tKe result 
of the virtual work statement (3.4) being written in terms of total stress and 
tractions rather than incremental quantities. Owing to this corrective nature 
and to the diminishing returns rne obtains from additional iterations, the 
second altematl/f. is used in the present study. 

Regul ation of Time Steps 

The creep calculations use a variable time step size which is auto- 
raatlcally regulated by the finite element program based on two criteria. The 
first criterion is the maximum percent difference between the incremental 


equivalent estimated creep strain and the incremental equivalent integrated 


creep strain for all the Gauss points in the mesh: 


= Max 




(3.17) 


The second criterion is the maximum ratio of incremental equivalent integrated 


This procedure could actually be more accurate if similar constitutive law 
integration procedures and equilibrium iterations are not performed with the 
tangent stiffness procedure. Also, it has been shown [4X] that many element 
types become overly stiff when using the tangent stiffness method for modeling 
constituc tve behavior approaching incompressibility. This problem is not en- 
countered with the current method. 
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creep atraiu to tho equivalent elastic strain: 



(3.18) 


Tho user specified, maximum pemtssiblc values for C and C are and C», 

i 2 1 ^ 

rospoctively. The size of the next step Is then obtained from 


At 


I+l 



(3.19) 


Note that the initial time increment cannot be determined from (3.19) and must 
be specified by the user so as to satisfy the two step size criteria. 

In the present study, the values of and are 0.2 and 1.0, res- 
pectively. With these values, it has been found that the initial time steps 
are controlled b’/ while later time steps are controlled by The values 

of and are strongly affected by the mesh refinement since a finer mesh 
results in Gauss points being closer to the crack-tip and therefore having larger 
stresses and strain rates. To determine the sensitivity of the solution to the 
selection of and C^, a compact specimen was analyzed with the above criteria 
and also with C. and C being halved (l.e., C •= 0.1 and C„ " 0.5), It was 

JU ^ 1 ^ 

found that the load point displacement differed by less than 0.5% for all 
time and that the steady-state solutions were essentially identical. It there- 
fore appears these values of and are small enough to ensure that the 
solutions to be discussed do not depend on these step size criteria. 
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SECTION IV 


ELEMENTS FOR SINGULAR CRACK-TIP BEHAVIOR 


Th.t8 section describes and compares several two-dimensional crack-tip 
singularity elements. Perhaps the primary motivation for introducing singu- 
lar crack-tip elements into finite clement models is the significant savings in 

computational expense. It is shown in Section V, for example, that 57 element 

- 1/2 

model with elastic, r , singular elements results in a more accurate solution 

than a non-singular 102 element model. The savings in GP time in this case is 

greater than 50%. Generally, one must consider that some additional effort 

is required to develop and implement a special crack-tip element and that this 

tends to offset the savings in CP time. It has been discovered, however, that 

the very commonly used, eigUt-noded, isoparametric element can be made to 
- 1/2 

produce an r strain singularity by merely shifting mid-side node locations 
via the node definition input data [42,43]. Therefore, a very convenient 
means for modeling linear elastic crack-tip behavior exists. It has also been 
shown [44] that a 1/r type strain singularity can be obtained with this element 
type thus providing a suitable element for non-hardening plasticity problems. 

For mote general singularity beliavlor, such as the r strain singularity 

associated with the HRR crack-tip field of power-law plasticity or creep, 
one resort to sepcially formulated elements. 

In the following sections, we consider special elements for linear 
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- 1/2 

eluatlc problem (r strain singularity) as well as special elements for 
problems involving the HRR, r strain singularity. In discussing 

these elements^ an attempt is made to point out their advantages and dis- 
advantages. 

Elements for Linear Elastic Materials 
Although many special elements have been used for linear elastic frac- 
ture analysis (see Atlurl (45) for review), we consider here only the eight- 
noded isoparametric element. There are two basic forms in which the eight- 
noded isoparametric element can be used as a linear elastic crack-tip element. 

In the first form, the two midside nodes adjacent to the corner node located at 
the crack-tip are shifted toward the comer node so that they occupy the quarter- 
point of their respective sides. This form is Illustrated in Fig. 4.1 by the 
Type A crack-tip mesh. The second form in which the element can be used is 
illustrated in Fig. 4.1 by the Type B crack-tip mesh, In this form, the elght- 
noded element is degenerated to a triangular element by defining two corner 
nodes and their midside node to be the same node which is located at the crack 
tip . Then the two midside nodes adjacent to the crack-tip corner node are 
shifted to their quarter-points. It is Important that only one node be used 
at the crack-tip, as opposed to three superposed nodes, since the latter case 
has been shown [44] to result in the 1/r type strain singularity. 

Barsoum [43,44] notes that numerical experimentation shows the degen- 
erate triangular form yields more accurate results than the nondegenerate 

element. He goes on to recommend that the four-sided configuration be abandoned 

- 1/2 

based on the premise that the r singularity exists only along the edges 

of the clement and also that the strain energy for this element becomes unbounded 

if exact integration is used [44]. Ying [46], on the other hand, concludes 
- 1/2 

that the r singularity does exist within the four-sided element as well as 
along its edges and that the strain energy for the element is bounded (and thus 
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Square 8 ~ node isoporatrietric Degenerate 8~node Isoparametric 
elements elements 


Fi.?. 4.1 Several crack-tip mesh configurations 



- 36 - 


the element stiffness is well defined). 

Apart from the above considerations, there are two aspects of the 
degenerate Type B element which inherently make it preferable to the nondegenerate 
element. The first is that the process of collapsing one side of the element 
to coincide with the crack-tip results in the element local coordinates being 
transformed into a form of polar coordinates. Since the element's shape 
functions are defined in the element local coordinates it is to be expected 
that angular bias will be much less apparent for this element type. 7he second 
feature of the degenerate triangular element which makes it preferable is that 
it is geometrically better suited for creating crack-tip finite element 
meshes of arbitrary refinement. Since the angular dependence of the near tip 
solution is significant, this flexibility for increasing the mesh refinement 
in the angular direction is important. Combining these two aspects with the 
numerical evidence cited by Barsoum, it seems the triangular, degenerate 
element is the better element for modeling the elastic crack-tip singularity. 

For this reason, all q; ssrter-point element calculations in the present study 
use the degenerate triangular form. In particular, the mesh configuration 
Type B of Fig. 4.1 has been used exclusively. 

Elements for Materials with HRR Crack-Tip Fields 

In the previous section, we discussed crack-tip elements for linear 

elastic material behavior. It has been seen that the standard eight-noded 

- 1/2 

isoparametric element can be made to have the r strain singularity and thus 
is useful for analysis of cracks in linear clastic materials. It can be shown 
(see for example, Atluri [45]) that this is the on]y singularity which this eight- 
noded element can exhibit. However, it can also be shown that higher order 
elements of the isoparametric family can result in singularities of the type 
^-n/(n+l)^ provided n is an integer. If we denote the order of the isoparametric 
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IntU'.rpo lac ion as m, then slnguIaL'lt*. ins of the type t)/t obtained for 

t an integer aiich that tsu [45]. In tenns of the HRR power law exponent n, 
this means we can have singularities of the type r where n in an integer 

such that n£m-l. 

Based on the above discussion, it can be seen that it is possible to 
employ isoparametlc elements as HRR crack-tip elements provided one is satisfied 
with integer values of the power law exponent, n. By choosing the highest value 
of n which ;>ne is interested in modeling, one can then program the n-fl order 
isoparametcic element. The problem with this approach is that values of n 
for common materials can be as high as 20. This implies that one would need to 
program an isoparametric element of order 21. While this is perhaps within 
reason, it will be slKU'm that nonisoparametrlc elements can be derived which are 
more readily implemented, 

Tw o Crack-Tip Elements from the. T.itcraturc 

In this discussion of special elements we limit consideration to two- 
dimensional, triangular elements with straight sides. The elements are derived 
in terms of the triangular polar coordinates (p,o) illustrated in Fig. 4.2 
and wliicli are related to the global cartesian coordinates 

“ xj' + p{xj - xj;) + |(a i- l)(xj - x^) } (4.1) 

In (4.1), the superscripts denote the node number. The crack-tip is assumed 
to be? located at node 1 (i.e., at p*0) . The geometric mapping of (4.1) is 
similar to circular polar coordinates in that the transformation cannot be 
inverted at ()“0. 

We now consider several choices Cor the assumed displacement fields 
within the triangular region. The first choice is 

u^(e,a) « u|^ + ^ (4.2) 
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IL can bo soon that (A. 2) Is Bliiiilar In form to (4.1) except that p Is replaced 
by so as to iiuhioo a sliigularlty in the displacement dertvativos . The 

throo-nodod elomont resulting from (4.1) and (4.2) is of the type proposed by 
Trucoy and Cook 147]. Inspection of (4.2) revouis that this element permits 
rigid body translation but does not permit rigid body rotations or constant 
strain modes, While problems exist for whleh this element provides reasonable 
resulLsi the lack of rigid body and o.onstant strain modes make this an un- 
desirable, elomont for general anlaysls. 

We now consider a straightforward procedure which allows an alternative 
to (4.2) to be written which (i) provides all the rigid body modes as well as 
a 1 1 the constant strain modes, (ii) results in the desired r'^ displacement be- 
liavior and (lii) results in compatible displacement fields with adjacent 
elements. First note that we can be assured our assumed dispiacemont field 
cunta.lns all rigid body and constant strain modes provided it can accomodate 
tlie following general displacement field 


"2 ■ “2 '’2’‘l ‘‘2*2 


(4.3b) 


where and are constants. Gleucly, a^ and provide for rigid 

translation modes, while b., and c,, provide for constant strains, t,, and 
respectively. The constants c^^ and b,, provide a rigid rotation if Cj^ « -b^ 4 0 
and a constant shear strain, •* b,, 4 0. If we substitute (4.1) 

into (4.3) and regroup terms we have the result 
Uj - a* + b*p + e*pO 


u,~ « a§ + b* 0 -f c$Po 

Ifc-N W ^ 


(4.4a) 

(4.4b) 


Starting from (4.4), we can now proceed to add terms us desired with the 
only condition being that we maintain compatible displacemonts with neighboring 
singular and nonsingular elements 
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We now cons Wcr a gentnvil appiroacli tot establishing displacemunt shape 
futuH lorivS for t riangular orack-tip eJt'mcnts. Since both displacement com- 
pum'nts will follow the same fornit we drop tlie nubscript fur simplicity. We 
now write 

u(.',u) « |(l + 0 )f.^(,v) (4.5) 


with 


r.,(i') » a, + b^P + c^p 


« a^ -f b^P + c.jp 


(4.6a) 

(4.6b) 


Inspection of (4.5) shows that on the element side l~2,u(p,-l) « f2(p) and 
on side l-3,u(p,l) “ also he seen that u(l,u) is linear on 

side 2-3. Since and each have three unlmowns, it follows that element 
sides 1-2 and 1-3 must have three. nodoi«. This means two new nodes must be 
created. Since the geometric properties of Che element do not depend on the 
locations of these nodes their positions along the edges of the element are 
arbitrary. In the following, however, we choose to place these nodes at the 
midsido.s. Tb.eso new nodes correspond to positions 4 and 5 in Fig. 4.2. De- 
noting the. nodal displacements by u'^ , j=l,5, we. now use the following con- 


ditions to determine the six unknowns in (4.6a) and (4.6b) 
u(0,o) « u'^; u(l,-l) “ u“; u(l,l) « u^ 
u(“-,-l) « i/*; u(i,l) “ u^ 

Tlio result is 


(4.7) 


«3 « uS b^ 


4 2 

2u - u 


2 1 , 3 

u - u - c„; b„ “ u 


1 

u - c. 


(4.8) 


.1-X 


* C 


„ 5 3 1 

2u - u - u 

2^*^ - I 


Substituting (4.8) and (4.6) into (4.5) and defining the functians multiplying 
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i ) 

u as we havo 


u(p,o) 


5 i 
X N^u 


(A. 9) 


where 


with 


.,5 T \ m5 i5,5 

Ni “ 1 - p - -jj(p - p); “ 'I’j'I'i; 


m5 .5,5 .,5 ,5,5 „5 ,5.5 

"3 *2'''2‘ "4 ■ V2’ "5 *2'>’2 


2l-‘ - 1 


t’ - f(l - a) 


>^2 “ + o) 


,5 

♦, - p 


g(i' - P) 


,5 2 , 4 . 

('2 - 8<P - P) 


By examination, it can be seen that (4.9) has terms similar to those of (4.4) 
and thus can represent al.l the rigid body and constant strain modes. From the 
form of (4.5,6j8) it can be seen that the element must be compatible with 
neighboring elements. Therefore, we have an element whihe satisfies all the 
roquircroents which we originally stipulated. If the a in Fig. 4.2 and in (4.9) 
is replaced by 2o' - 1 we recover the form of the equations suggested by Stern 
(481. 

As a result of the appearance of in the shape functions for this 
element, the integration involved in evaluating the stiffness matrix (particularly 
the integration with respect to P) is not suited to Gauss quadrature. In [48], 
Stern derives a special Integration rule, which when combined with standard 
Gauss quadrature for integration with respect to a, results in exact stiffness 
integrations. Unfortunately, the expression (28) in [48] which specifies 
the relative radial location of the quadrature points is apparently in error. 

The corrected expression is as follows: 
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^.L “’'2 ^ X)4"-^J 


(4.10) 


where Xj anr^ x,, Uenute the values of |i at which the sampling points are located. 
The corresponding weights are then given by 


1 X-l 

2 ^2 


1 

X+1 


\-l X 
^’^2 - ’‘l 


w 


1 X-l 

'2 ^^2 


\+l. 


2 X~1 X 


(4.11) 


(4.12) 


^ 2 X— 1 

This quadrature rule integrates tenns of the type t"* , 0 and p ~ ‘ exactly. 

Since the rule has four parameters (x^,x^,w.j ,w^) and is only required to integ- 
rate tlu'ce types of terras, the locations of the quadrature points are not 
uniquely defined. Selection of x.^ according to the following cirtcrion resulta 
in both the numerator and denominator of the bracketed terra of (4.10) being 
positive and thus results in a valid quadrature rule. 



(4.13) 


Stern (48] presents a family of elements which are developed so as to 
bo compatible with surrounding nonsingular polynomial based elements of arbit- 
rary order. While it should be possible to verify that each member of this 
family dues indeed satisfy the requirements which were discussed in deriving 
the above element, this procedure provides little Insight to Che method for 
deriving such elements. In the next section, we generalize the procedure used 
in arriving at (4.9) to derive an element which is compatible with quadratic 
eleraunts (e.g., eight-noded isoparametric elements)* 

Derivation of a New Crack-Tip Element 

In this section we generalize the procedure used in the previous section 
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to dorlvo a crack-tip elumcMU which (i) contains all rigid body and constant 

strain modes, (ii) results in compatible displacomcmt fields with neighboring 

singular and nonsingular elements and (iii) results In an atbitrary strain singu- 

\-l 

larity at the crack-tip fo the type r 

We start in tlve same manner as before by writing the displacement field 
in the form 

u(p,o) » (l-a^)E^(p) + (4.14) 


who ire 


f,(.» 


"i + '■i” 


+ c.n 
i 


(4.15) 


Noting that the form of (4.15) requires three displacement*: being specified 

along each of tliree radial line segments (o** -1, o *=> 0, o “ 1) we introduce 

node points at locations 4 througli 7 as iriustrated in Fig. 4.2. At this 

point it is seen that this procedure will result in an interior node. It 

will be shown later that this node can he eliminated in a number of ways. 

i 

Denoting the nodal displacements by u' , J“l,7 we now use the. following 
conditions to detennlne the nine unknowns in (4.14,15): 


u(0,o) “ u^; u(l,-l) ** ; u(l,l) = u^ 

u(“-,-l) “ u^; u(’|,l) >= u(l,0) >= u^ 

i»(Y, 0) “ u^ 

The result is 

1 

* ^2 ” ^3 “ “ 

, 6 1, 2 1 ,31 

*» u - u - c.^; b, » u - u ~ c^; b^ “ u - u - 

i/o 7 6 1. J.._ 4 2 1. 

^1 ” “ *■' - ) 5 - u - u ) 

i/o 531. 

C.J - -(2u - u - u ) 


(4.16) 


(4.i7) 
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, > > .,1-^ 
wlLh « 2 - 1 

Substituting (4.17) and (4.15) Into (4.14) and defining the functions multiplying 
u as N, we have 


u(f),u) 


7 

£ 

i-1 


J i 

N u 


(4.18) 


whore 


nJ. 

1 - p 

1 

1 

7 

.7 , 7 

7 77 



N 3 “ '<*3'J'2‘ 

7 

.7,7 

7 7 7 

n; « 
5 




♦ 2*1 


. 7,7 


w it li 


=> 1 - 0^ ; ,{'2 “ “ -^) 5 'J’3 ^ 


,7 2, A , 

i|'^ « -jCp - p) 


1>2 " |t(l + 3)p - P^J 


Inspection of the derived shape functions verifies that this element satisfies 
all the requirements which we stipulated at the beginning of the derivation. It 
was noted that this element has an interior node. While interior nodes are 
generally avoided so as to reduce the bandwidth of the equations to be solved, 
it seems that in the case of crack-tip elements the advantage of having additional 
degrees of freedom in the vicinity of the crack-tip more than compensates for 
the few additional equations which are involved. 

We now consider several alternatives for eliminating the interior node 
of tills element and note that one of these results In the corresponding element 
of Stem [48). We start by substituting (4.15) into (4.14); using the condition 
that * f 2 ^^) “ f.j(0) and regrouping terms we can write 

u(p,a) » {a + b^p + Cj^p^) (con't on next page) (4.19) 
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It should be uadefstood that any arbitrary condition relating u to one, several 
or all of u through u will suffice to remove the interior node (i.e., node 7) 
from the clement relations. However, it seems more natural to eliminate the 
node by the removal of one of the terms of (4.19). Noting that the first 
two terms of the first square bracket and the first term of the second square 
bracket represent the rigid body modes and constant strain modes, we are left 
with four terms which can possibly be deleted so as to eliminate the interior 
node. 

If we choose to eliminate the term hy constraining Cj^ to be zero 

wo see from (4.17) that this implies u^ = ^ choice would some- 

what defeat Che objective of having singular displacement derivatives- and thus 
is not advisable. Furthermore, it is inconsistent to retain the higher order 

terms op'^ and while not retaining p^. Note that the term op^ cannot be 

7 7 

used to eliminate u since u does not appear in its coefficient (i.e., in 

2 2 V 

either c^ or c^) . Therefore we are left with the terms o p and a p ^ Either 

7 

of these terms can be chosen to eliminate u , Stern's element [48] corres- 

2 . 

ponds to the case in which the coefficient of a p is set identicaxly to zero. 

Of the elements discussed above, only the seven-noded element has been 
implemented in the present study. In Section V, this element is used for 
elastic analysis as well as for creep analysis. The special quadrature rule 
proposed by Stem [48] and summarized in (4.10) through (4.13) has been used - 
exclusively in evaluating this element's stiffness. 




SECTION V 


CREEP CRACK GROWTH COMPUTATIONS 


Description of Problems 

Tl>e creep crack growth analyses which will be presented in this chapter 
deal with three distinct problems. The following sections introduce each 
pfjblem by describing the physical aspects such as geometry, loading 
and material properties as well as by describing why the problem was selected 
and what is hoped to be gained by its consideration. 

All calculations in this chapter assume Infinitesimal strains and small 
deformations. The crack propagation calculations use quarter-point crack-tip 
elements and a mesh shif ting/remeshlng procedure. 
problem I; Non-Steady Creep of a Compact Specimen 

The compact specimen geometry was chosen for study because of its wide- 
spread use in fracture experiments and because numerical solutions for this 
problem have appeared in the literature thus providing results with which to 
compare. The dimensions of the specimen as well as the material properties 
and applied loading (see Fig. 5.1) were chosen to coincide with those used 
recently by Ehlers and Riedel [33]. The problem is used for a mesh refinement 
sensitivity study and for exploring various aspects of the con- 

tour integrals during both nonsteady and steady-state creep. 

Several finite element meshes have been uased in the analysis. All 

- A7 - 




oi‘ thcsu meahCB employ two-dimensional, elght-noded, Isoporaroctrlc elements. 

The integrations for these elements arc accomplished with 2x2 Gauss quadrature 
and therefore only elements with straight sides are employed. As seen from 
the meshes in Fig. 5 , 2, the pin-loading hole is not modeled. In all models 
the horizontal placement of the point load corresponds with the load line of the 
ASTM standard geometry (x ■ 25.0 mm). The vertical position is y - 32.5mm. 

A sensitivity study showed that shifting the load to y ■ 40 mm has virtually 
no effect on the pertinent aspects of the solution. 

Most of the meshes contain collapsed eight-noded isoparametric elements 

at the crack-tip as illustrated in Fig. 6.2. In several calculations, the 

midside nodes of these crack-tip elements are shifted to their quarter-points 

- 1/2 

so as to produce an r strain singularity at the crack-tip. Also, several 
calculations are performed with a special conforming seven-noded, triangular 
element which Imposes the HRR, r type strain singularity. Table 5.1 Iden- 

tlCifo the meshes for which calculations are made and also gives the load point 
displacement and for the elastic solutions. These values are compared to 
those based on the expression for given by Srav^ley [49]. 

Problem II; Constant Velocity Propagation in a Creeping Strip 

This problem is concerned with a finite height, infinitely wide strip, 
with a semi- inf Inite crack. Loading consists of uniformly applied displacement 
rates at the top and bottom edges. This problem has been chosen for two 
reasons. First, since the strip is Infinitely wide and the boundary conditions 
do not change with time, the propagating crack-tip fields can be expected to 
reach a 'convecting steady-state'* creep condition. Here we use the phrase 
"convectlng steady-state" to mean that the field remains unchanged in time with 
respect to a coordinate system which is centered at and moving with the crack- 
tip. This terminology is used so as not to confuse this condition with the usual 
steady-state creep condition in which material stress rates are zero. 
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55 elements’ 


OF POOH hH* ‘ ‘ 


198 nodes 
384 d.o.f, 

57 elements i 

200 nodes 
388 d,o.f. 













IT 

♦ 








1 

J 


t 

"“"1 








f 

1 

1 

4 

I 







t 



) 

1 

1 

r - 
1 

— 


— 

- - • 

» 

1 

t 

t 

1 . 

1 

1 

1 



i 

t 

- --L. , 

t 

1 

1 

r- 

i 



1 

1 

t 

t 

1 

1 

1 

1 



1 

1 

i 

I 

1 

1 

1 

1 

s 

\ 

/ 

/ 

1 

1 

« 

1 

* 

1 

— u. 

♦ 

1 


0. The 55 (square crack tip elements) and 57 (triongular 

crack tip elements) element mer' es 



b. The 102 element mesh (331 nodes; 642 d,a.f.) 



Fig. 5.2 Finite element meshes for the compact specimen 

(contour integral paths are indicated by dashed lines) 
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Table 5.1 Suininary of Computational Aspects and Comparison with Results from the Literature 
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stopped at 600 hours; solutions 9-10 arc stopped at 300 


In tho convecCing stecdy-state case, stress rates for material points 
are not zero. .\s a result, C* cannot (In a strict sense) be a valid crack-tip 
parameter. As should be clear from earlier discussions, • valid 

parameter at all crack speeds since It does not require material stress rates 
to be zero. From a practical point of view. If the crack speed is low enough, 
then one can expect C* to be a useful quantity. By varying the crack propagation 
speeds for this problen over the range of velocities observed experimentally 
(for a given material). It should be possible to determine If this range has 
any portions In common with the range of velocities for which Is a useful 
parameter. 

The second reason for choosing this problem is that C* can be evaluated 
aiialytically for the special case when the crack Is stationary. This allows an 
independent check on the finite element calculations and serves as a reference 
for the analyses In which the crack is propagating. The analytical evaluation 
of C* parallels the evaluation of for a similar elastic strip problem as 
discussed by Rice [50]. (See Appendix E) It should be noted that C* has been 

a 

shown to be related to the steady-state value of (T.) and therefore it is 

1 c 

possible to obtain stationary crack case from C* and equation 

(B.l) of Appendix B. The direct evaluation of in terms of either Its in- 

tegral representation or Its energy representation requires knowledge of the 
stresses in the region of the strip adjacent to the crack-tip and therefore Is 
not a trivial task. 

The material properties used in this problem are representative of 30A 
stainless steel at 630^C. These material properties and the finite element 
discretization are given In Fig. 5.3. The mesh for this problem may at first 
appear rather coarse; however, elastic and steady-state creep solutions obtained 
with this mesh are sufficiently accurate to Justify its use. The comparison 
of computed elastic values and steady-state C* values with their analytic 
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Fig. 5.3 Strip geometry , finite element mesh and material properties 
(contour integral paths are indicated by dashed lines) 
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valuos is giViMi in Tiil)U' 5.2. 

i i i > C»‘ o o.p Crack Gi't wth in Doub ln-MdRu>-Cy.nol- iiJlO£.iiiu'im 
Tin? inirpoHO of ('on,sitlo.Tl.n8 tills pfoblom is l.o upply the modol to a 
pt’obU'm fui' wliicti oKpcrimontal data exists. Wlvtie mueli experiment. a i data has 
been reported in the, literature, most autliors do not include sufficient in- 
formation to allow a numerical simulation of their experiments. Tlic current 
problem is based on tin? exporimcnts of Koteraxawa and Twata [51]. Tlie primary 
reasons for selecting this work for study are that crack length versus time 
liistories were given and that the experiments were porformed on 304 stainless 
steel foe wiilch high temperature elastic iUid creep properties wore already 
available. 

The 8eometr\r of the experimental speciinens is given in Fig, 5.4, The. 
finite element mesh for the ealcul,itions i.s shown in Fig. 5,5 with contour 
integral paths being Indicated by dashed lines. It can be seen that the. mesh 
takes advantage of tlve two planes of symmetry for the specimen and does not 
model the 60%oteh.^‘ The initial crack length indicated in Fig. 5.5 corresponds 
to the. notch depth in the, specimen. All calculations for this specimen assume 
plane stress conditions and use tlm material properties given in Fig. 5.3. 
Elastic J,^ results for two crack lengths are compared in Table 5,3 with those 
(based on formulas for Kj) from [52] and are seen to be in good agreement. 

Compact Specimen Analyses 

The following describes several calculations for a compact specimen 

during transition from .m initial elastic state to one of steady-state creep. 

The geometry j loading, material properties and other details were described in 

the first section of this chapter as Probelm I, We first consider results for 

* 

the 300 clement mesh of Fig, 5.2 in terms of Ct and then (T,) . Then we address 

.. , . X 1 c 

Hlodeling the notch would have required the. mesh shifting subroutines to be 
goneralised. 
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double -edge -notch experiments 
(contour integral paths are indicated by dashed lines) 
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Table 5.3 Computational Aspects of the Elastic and Non -Steady Creep 











tlie toptua of nu'sh reflnerount and the use of spectel crack-tip singularity ole- 


men ta. 


Ca Lculati.on of Nonsteady Creep 

Thie. path-indepondance of during nonsteady creep is illustrated 

in Fig. 5.6 using results from the 300 element tnesh. The f, superscript 
designates the particular ^234 contour which is used, with being the 
non dimensional distance from the crack-tip to the point where the contour 
crosses the erack plane. Therefore, fj is zero at the crack-tip and has a 
maximum value of unity when the contour is at the boundary of the specimebi. 
Values of (C|) ' are plotted as a function of time for nine values of S ranging 
from 0,03 to 0.92. It is seen that ' is largest for contours close to the 
crack-tip (small ?;) and that as steady-state is approached, the values from 
ail contours converge to C|. The solution has essentially reached steady- 
state at 300 hours. After 300 hours, the values of (C*)^ for all nine contours 
arc within 1.5 percent of their average value. This value of C^, as well as 
values from calculations with the other meshes, is given in Table 5.1. 

Now we consider computed values of (T-) as approximated by (AT.) /At. 

jL c r c 

The values of (AT are obtained through the specialization of (2.10) to the 
X c 

case of infinitesimal strains, small deformations, symmetric mode I be' avlor 
and traction-free surfaces: 
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Since (T,)^ is the limit of (T^) as c goes to zero, (T.)^ is plotted as a 
function of e for several times (sec Fig, 5.7). In this figure, e is the 

nondimensional size of V and is measured in the same manner as X, the 
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time (hours) 


(C^)^ as a function of time for 
1 contour integral paths 
m 300 element mesh) 









nondimens I on a I size of T, 


The open points ace the values of eom- 


putod by (5fl), fot nine eontours in the 300 element model. The value of the 
ccack-tip parameter (T.) is given by the intersection of each respective curve 
with the t; * 0 axis. Due to the large gradient in for small c it is 

seen that the accuracy of any extrapolation based solely on these evaluations 
of (T*) (i.c., open points) would be of questionable accuracy, except perhaps 

4- C 

near steady-state conditions. The solid points at c ■ 0 in Fig. 5,7 have been 

obtained using (2.14). It is seen that these values of (T.) appear to be reasonable 

X o 

extrapolations of the curves of (Tj^) (5.1) thus giving some degree of conf idence in their 

accuracy. 

Path Independence of (Tj^)^ 

Based on arguments put forth in earlier portions of this paper, the 
value of obtained through (2.14) should be independent of the path 

(i.c. , Which is used in its computation. This path-independence is 

illustrated by Fig. 5.8a where (Tj^)^ is plotted as a function of the nondimen- 
siunal distance of ^234 from the crack-tip, for several times. Generally, 
the path-independence is seen to be quite good. The largest deviation from 
path-independence in this figure is for the Intermediate time of 10,8 hours 
with the difference between the extreme contour values being less than three 
percent. To further emphasize this path-independence, (Tj^)^, Is plotted as a 
function of time in Fig, 5.6. As a result of its path-independence, (Tj^)^ is 
represented by a single curve. Interestingly, this curve is a straight line 
for times before approximately 10 hours. 

Riedel and Rice [36] have arrived at the following approximation for 
K (which they call A(t) based on the assumed approximate path-independence of 


during the initial portion of nonsteady creep: 
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Comparing (5.2) with (2.50) one concludes that (tj^)^ should behave like 1/t for 
times wlten (5.2) is valid. In a log-log plot of (Tj)^ versus time this would 
result in a straight line with a slope of -1. The straight line shown in Fig. 

5.6 Is inclined from the horizontal by 40° and therefore has a slope of -0.84, 

The current work has resulted in some evidence that is approximately path- 
independent during Initial nonsteady creep but that its value tends to increase 
with time. This tendency for to Increase with time could explain the rather 
significant departure of the current results from the behavior of (5.2). 
Quarter-Point Singularity Element Calculation 

We next consider the results of computations using 57 and 102 element 
meshes with quarter-point singularities. The purpose of considering these less 
refined meshes is to determine if the expense and effort in using the 300 element 
model is necessary for obtaining accurate results. Table 5.1 summarizes the 
results of these meshes for the limiting cases of purely elastic behavior 
and steady-state creep behavior. For the elastic problem it is seen that the 
results from these meshes agree with the 300 element mesh results to withi'^: 
one percent. At steady-state the 102 element model still agrees with the 300 
element mesh (in terms of C^) to within one percent while the 57 element 
model now differs by approximately eight percent. 

The contours used for the 57 and 102 element mesh are indicated in Fig. 
5.2. The 57 element mesh has four contours while the 102 element mesh has 
eight. The path-independence of (fj^)^, as computed frok '2.14), is illustrated 
for these two meshes in Fig. 5.8b and 5.8c. It is seen that the degree of 
path-independence in both is similar to that observed for the 300 element mesh. 
Since we have evidence that Che 57 element mesh is less accurate than the other 
meshes at steady-state, it appears that high quality of the path-independence 
cannot generally be interpreted as meaning the solution is accurate. 



To determine the adequacy of the 57 and 102 element meshes for the non- 
steady creep problem we now compare their histories with that obtained 

with the 300 element mesh (see Fig. 5.9). The curve appearing in this figure 
has been placed through computed points from the 300 element mesh. The results 
of the 102 element mesh agree almost perfectly with this curve for times between 
0.2 hours and 16 hours. Prior to this period and after this period the results 
fall below the curve by as much as 20 percent. While little can be said about . 
the absolute accuracy of the calculations for early portions of nonsteady creep, 
we know (based on Appendix B) that (T^)^ should agree numerically with C* at 
stejidy-state to within a few percent. Therefore it -can be said that the 

values of (T.) from the 102 element mesh are significantly in error at steady- 

1 c 

state. Recalling that this model gave a steady-state value of C* which ®gree 
quite well with the 300 element mesh results (see Table 5.1) it is perhaps sur- 
prising that such a significant error in the steady-state value of -Can 

exist. To better understand the results of this model, is plotted as a 

function of ^ in Fig. 5.10. It is noted from this figure that the values of 
(T^)^ based on (2.14) (i.e., the solid points) appear to be reasonable extrapo- 
lations for tiroes when the results are in agreement with the 300 element xnesh 
results. However, as steady-state is approached, it is seen that these solid 
points no longer appear reasonable. If, however, one extrapolates the values 

of (t-)^ to e = 0 for the bottom two curves of Fig. 5.10, it is found that 

1 c 

these values of (T.,) are in good agreement with the 300 element mesh results. 

1 c 

In comparing the equations for evaluating C^, and it is 

. 

seen that (T.) is the only one of the three which involves an integration 

* X c 

over the crack-tip quarter-point elements. Based on this and the apparently 
good accuracy of and (ij^)g It is believed that the integration over 
these elements is the major cause of discrepancy in (Tj^)^ between the 102 
and 300 element mesh calculations. 
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Fig. 5.9 Comparison of Ci]_)Q from caloulations 
witjh three finite element meshes 
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Tl\e 57 element iresulta do not compare favorably with the curve of Klg. 

5.9 for any significant portion of the solution. For most times the values of 
fall below the. curve with the percent difference ranging from 50 percent 
at 'c •* 0.02 hr. to 1.5 percent at steady-state. Based on the discrepancy of Cj^ 
indicated In TabXii 5.1 and in the generally bad comparison of In Fig. 

5.9. it appears that the 57 element mesh with quarter-point singularity is not 
sufficiently refined for accurate creep calGulations. thl.s conclusion is per- 
haps a bit unexpected considering the degree of accuracy which this mesh displayed 
for the elastic probelm (see Table 5.1). The reason for this drastic change 
of accuracy in going from elastic to creep behavior may be that the crack-tip 
strain singularity (i.e. , is inappropriate for the type be- 

havior expected to exist during creep. This topic is addressed in the following. 
HRR Sin gularit y Element Calculations 

Based on the above observations, several analyses have been made using the 
sevon-noded variable singularity clement described previously (Section IV). The 
clastic solutions obtained with this element agree very well with those using 
the quarter-point isoparametric element as can be seen from the entries of 
in Table 5.1. Also included in Table .5.1 arc t!»e quasi-steady-state values 
of C^. The introduction of the correct strain singularity for steady-state 
creep (r does not significantly affect the 102 (Element mode's CJ but 

does Imporve that of the 57 element model. 

The analyses which use the seven-node singular element have the same 
singularity for the elastic so.lution and the subsequent creep solutions. At- 
tempts at changing the singularity from the clastic to the 

value between the elastic and first creep solution have created numerical dif- 
ficulties due to the disequilibrium introduced in the process. No attempt 
at a gradual transition has been made. 
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The results using the. seveir-noded slngulnr oleinent are shown In 

Fig. 5. 11. The solid curve represents the results of the 300 element model. 

The evalvuition of according to (2.14) with the numerical procedures 

being identical to those employed with the quarter-point elements except for 
the contribution of the singular elements to the vo.lume integral. For the 
quarter-point elements* the stresses are assumed to bo distributed lineraly with 
respect to the local coordinates. The volume integral is then evaluated in 
terms of quantities at the 2x2 Gauss ponts. For the seven-node element, linear 
interpolation is used and in addition, several calculations are done assuming 
radial dependence of the type o *> a + b . It can be sosn from Fig. 5.11 

tl\at none of the calculations agree well with the 300 element results. 


Based on this set of calculations, the general disagreement In (T.) 
between the singular crack-tip element models and the 300 element model does not 
appear to be due to the strength of the singularity which is introduced at 
the crack-tip. The general accuracy of G|’ for all the solutions with either 
elastic or creep type strain singularities supports this view. Rather, it seems 
likely that the difficulty in computing the volume integral of (2.14) stems from 
the crack-tip element fields not satisfying the condition 


Lt 

e-k) 



(5.3) 


From the discussion of Appendix A it can be seen that if this condition is not 
satisfied while at the same time the fields have the correct asymptotic (singu- 
lar) radial dependence, then the volume Integral of (2.14) does not exist. 

It therefore appears that accurate evaluation of ('f^)j, asing (2,14) 
cannot be accomplished if one uses crack-tip singular elements which provide 







for the satisfaction of condition (5.3). At this point, it appears that one must 
use a rather refined non-singular mesh (such as the current 300 element mesh) 
or Introduce special crack-tip elements which satisfy (5.3) in order to compute 

(T.) accurately. The next section Illustrates that for many probelms of 

♦i. c 

practical importance, a more attractive alternative may exist, 
as a Crack-Tip Field Parameter 

The previous discussion has pointed out some computational difficulties 

involved with evaluating (T ) . It was concluded that these difficulties 

1 c 

are associated with the contirbutlon of the crack-tip singularity elements to the 

volume integral of (5.1). It has been seen that despite the elastic strain 

singularity Introduced by the quarter-point element scheme, the 102 element 

mesh gives accurate values of C*. Assuming this reflects the general accuracy 

of this solution, it is desirable to use this relatively inexpensive model as 

opposed to using a very refined non-singular mesh or to introducing a special 

crack-tip element which satisfies condition (5.3). 

The effect of deleting the crack-tip singular elements from the volwrae 

integral of (2.1A) is shown in Fig. 5.12. Deleting these elements means that 

we are In fact evaluating (T.) where V is the volume encompassed by the 

X c e 

* e * cS 

crack-tip elements. We will denote this particular (Tj^)^ as (Tj^)^* It will 

be shown that depending on the relative size of the crack-tip elements and the 

• 6 

proximity of the solution to steady-state condition, (Tj^)^ is a good approxi- 

mation to (T-) . 

1 c 

The solid curve in Figv 5.12 represents the results of the 300 element 

* *5 

mesh. The dashed curves are (T,) in the case of the 57 and 102 element 

1 c 

meshes and is (T,)^ with e = 0.03 in the case of the 300 element mesh. The 
J- c 

crack-tip element sizes for the 57 and 102 element meshes are 10 mm and 2.5 mm 
(or 20 and five percent of the ligament size), respectively. 
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3CX) (t,)c 

300 percent) 



Fig, 5.12 Comparison of and histories 

for several finite element meshes 
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6 

The values of (T. ) for the 102 element mesh cotneide with the solid 

X C 

cuve for times after about 30 hours* Therefore (T^) Is a velldj path indepen-> 

JL O 

dent, crack-tip parameter for times after about 30 hours and for values of (T*) 

C 

beginning at approximately 1.6 of the steady-state value. Fig. 5.6 shows that 
is still significantly path- dependent at 30 hours and thus is not an acceptable 
crack-tip parameter until much later. 

The curve of (T,) , (e ■ 0.03), for the 300 element model seems to Indl- 
1 c 

cate that the validity of (T.) can be expanded to earlier times by reducing 

JL C 

the size of the quatter-po‘'nt elements. For example, a 6 of three percent of the 

. g 

ligament would apparently result in being valid as early as seven hours and 

for values of (T.) as large as 4.3 the steady-state value of (T ) . The curve of 
X C X C- 

(T^) ^ for the 57 element mesh tends to approach the solid curve as steady-state is 

approached but never acutally converges even at steady-state. This indicates that 

• ^ 

this mesh is too coarse for fn be a useful parameter. 

Constant Velocity Propagation in a Strip 
We now present some calculations for the cracked strip problem previously 
referred to as Problem XI. The geometry, loading and material properties for 
this problem are summarized in Fig. 5.3. The purpose of this problem is to 
determine how significantly the crack-tip fields are affected by crack propa- 
gation velocities commonly observed 1.) experiments. If for realistic crack 
speeds, the crack-tip field is essentially the same as for a stationary crack, 
then C* is path Independent and characterizes the crack-tip fields. In any 

case, (1’ ) is a valid parameter. 

X c 

As noted previously, the steady-state C* values for the infinite strip 
problem can be obtained analytically without much difficulty. (See Appendix E.) 
Therefore, the procedure for this set of calculations is to select three values 
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of C* which span the* range of values reported in the literature. The values 
which have been chosen are 0.05, 5.0 and 50 N/mm,hr. Having these values, the 
corresponding remote steady-state is determined as well as the edge dis- 
placement which will result in the same remote elastic These displacements 

are applied to the model elastically at t - 0. The resulting values of are 

compared to the analytic values in Table 5.2. Next, the steady-state edge dis- 
placement rates are determined. Using the elastic solution as an initial 
state, the displacement rate, 5, is applied until the model reaches steady- 
state, The computed steady-state values of CJ are compared to their analytic 
values in Table 3.2. The next step is to determine an upper bound crack 
velocity for each of the chosen values of CJ. The following formula is based 
on the experimental data reported in [23,241 and represents data from center- 
crack, double-edge-crack, single-edge-crack, compact, and round-bar specimen 
types . 

^ ^ oCJ (5.4) 

where 

_2 

1,68 . 10 (upper bound) 

a = 

3.36 , 10 (average 

Having reached steady-state, the crack Is propagated at the upper bound velocity 
given by (5.4) until it is determined that a convecting steady-state has been 
reached, 

As noted previously, these calculations use the quarter-point crack-tip 
element. The crack growth simulation is accomplished through a combination of 
mesh shifting and remeshing as described in Appendix D. The nominal size of ' 
the crack growth Increments is 0.4 mm or two percent of the crack-tip element 
width. For the highest velocity case (Gj ■ 50 N/mm.hr), this procedure re- 
sults in crack growth at approximately every fifth solution step. 
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Reaulta for a Plane Strain Strip 


Tl»e results of the plane strain strip calculation with Cj ■ 50 N/nun,hr 

and ^ * 1.65 mm/hr are given in Fig, 5.13, The values of (Tj^)^ and C* are 

given for the portion of the calculation prior to steady-state as well as 

during the crack propagation protion. The band represents the range of values 

obtained from the four contours illustrated in Fig. 5,3, Both and CJ 

converge to the 50 N/mm.hr value at steady-state. During the crack propagation, 

it is seen that (T.) and G* do not depart significantly from their steady- 
«l« C X 

state value. This means that this combination of loading and crack spr^ed 

results in the crack-tip fields being essentially at steady-state conditions. 

This in turn means that both (T,) (or (T,)„) and C* are valid crack-tip 

xc xc X 

field parameters. 

A closer view of the crack propagation portion of these curves is given 
in Fig. 5.14. The dashed curves bracettng the initial portion of the solid 
curves represent the degree of path-independence and continue to be representative 
of the path-independence observed during the crack propagation steps. For both 

(T.) and C*, it is seen that the strip has essentially returned to Its steady- 

1 c 1 

state condition prior to each crack growth increment. It is tViought that the large 
departure of (T. ) (as compared to Ci*<) is more representative of the nonsteadiness of 

X O 

the crack-tip field, since the validity of C| in general, and the numerical evalua- 
tion of W*** (2.43) in particular, are based on the existence of steady-state conditions. 
The effect of remeshlng is seen at approximately eight hours. The 
first two steps after the remeshlng were found to result in rather erratic con- 
tour integral values and are not indicated in these figures. The equilibrium 
correction feature of the present model and the automatic time step regulation 
procedure both act to quickly restore equilibrium at the crack-tip. 

The propagation portion of the calculation with CJ ■ 5 N/mm,hr and 
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da 

dt 


» 0.1X1 mm/hr Is given in X'ig. 5,15. Her© again it Is seen that both 

and hav , converged to the analyticaX value of (to within two percent, which 

Is also about the degree of path-dependence) . Comparing these results with 

those in Fig. 5,14 for the higher C* and crack speed It Is seen that steady- 

state creep conditions were not reached until 12 hours as opposed to approximately 

two hours in the previous case. Also, the return to the Bteady-state value after 

mesh shifting takes more time (two hours compared to 0.25 hours). However, when 

compared to the time between crack growth steps (both use, 0.4 mm) it is seen 

that the lower velocity case return to steady-state well before the next growth 

ste,p occurs. This result indicates that lower load levels and crack speeds 

are Inherently closer to steady-state conditions. While this behavior may 

seem intuitively correct, it should be kept in mind that these results depend 

on the empirical formula (5.4) which is only valid for 304 stainless steel. 

It remains to be seen if simf.lar behavior occurs in other materials. 

A calculation has also been done for the case of C^“0,05 N/iiun.hr. As 

a result of the large number of solution steps between crack growth steps, when 

using the maximum velocity of 5 . 10 ^ mra/hr, the calculations used a higher 
-3 

velocity (5 . 10 mro/hr) . Even at this unreallsltcally high velocity (for this 

level of loading), the behavior is more steady-state-llke than the case of C|»5.0N/mm. 
hr described above. 

Results for a Plane Stress Strip 

In both plane strain problems discussed above, the steady-state value of 

(T,.) is equal to C* to within the accuracy of the calculations. This is con- 
X ^ X 

sistent with the relationship aivd comparison of C* and steady-state (T. )^ given 

X X ’Ci 

in Appendix B. According to the approximate numerical values of this appendix, 

* 

there should not be as close agfeetnent between C| and (Tji)j. iu the case of 
plane stress. The primary purpose of this analysis is to verify this predicted 
behavior. 
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(C| = 5.0 N/mm • hr ; v * o.lll mm/hr) 
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For this plane stress analysis, C* was chosen to be 50 N/mm . hr and 
the crack was again propagated at 1.65 nun/hr. The remote t , the steady-state 

yy 

displacement rate, and the elastic displacement, 6, are 171 MPa, 0.168 mm/hr 
and 0.114 nun, reapectlvely , 

The results of this calculation are given in Fig. 5.16 and 5.17. It 

.6 

ib seen from these figures that (T, ) does converge to a somewhat higher 

i c 

value at steady-state then CJ. The steady-state Is seep from Fig. 5.17 to be 
approximately 52 N/mm . hr which is higher than Cj by four percent. While 
this is a somewhat smaller difference than suggested by Appendix B, the sign of 
the difference is the same. In light of the approximate integration used in 
obtaining the numeric values in the appendix, this discrepancy is within reason. 
As expected, the general behavior for plane stress conditions is essentially the 
same as for the previous plane strain analyses. Therefore, previous observations 
concerning the steady-state nature of the crack-tip field during crack propa- 
gation are unchanged by the shift to plane stress conditions. 

Double-Edge-Crack Specimen Analysis 

The following describes *<everal calculations and their results for the 

probIcM previously referred to as Problem III. The geometry and finite element 

mesh for the double-edge-crack specimen are given in Figs. 5.4 and 5.5, respec- 

o 

tively. The material properties are those of 304 stainless steel at 650 C and 
are assumed to be the same as those used in the strip analyses. (See Fig. 

5.3) Calculations have been made for remote applied stresses of 157 and 176 MPa. 
The experimental crack growth histories for these two stress levels are repro- 
duced from [51] in Fig. 6.18. It Is seen from these curves that the first two- 
thirds of the specimen lives are characterized by crack velocities of less than 
0,01 iran/hr compared to nearly 0.5 nun/hr as rupture is approached. 

The primary purpose of the following calculations is to verify the 
conclusions which were reached in the previously described strip calculations; 
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Fig. 5.17 crack propagation portion 

of a plane stress strip analysis 
(CJ = 50 N/mm • hr v = I.65 mm/hr) 
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Fig. 5.18 Double- edge- crack specimen crack growth histories (from [773) 





that is, that th«i crack-tip fields are essentially creep steady-state fields 
even for the most rapid creep crack velocities. These calculations will be a 
valid check because the input to the calculations is only the remote applied stress 
and the measured crack velocity history, and does not in any way depend on experl- 
mental determination of C* or ‘lid strip calculations. In fact, 

Koterazawa and Iwata do not report such measurements in [51). 

Analysis of Initial. Low Velocity Crack Growth 

This section describes the simulation of the initial portion of the 
crack velocity histories given in Fig. 5,18. In all of these calculations, the 
entire load is applied elastically at t “ 0 and held, constant throughout the 
subsequent creep solution steps. The convergence of and C* to their 

steady-state values is shown in Fig. 5.19, with the dashed lines in the C* 
plots denoting the degree of path-dependence. It is seen that steady-state 
conditions are reached between a half and one hour after the load is applied. 

(Table 5.3 summarizes the computational aspects of this portion of the calcu- 
lation.) Therefore, it is seen by refereing to Fig. 5,18 that crack growth does 
not begin in the two specimens until well after steady-state conditions are 
reached. Since the current calculations assuiae small displacements and infinitesi- 
mal strains, and since only the strain and displacement magnitudes depend on time 
once steady-state is reached, there is no reason to continue the numerical calcu- 
lations to the crack initiation times indicated by the experimental results. 
Therefore, the initial crack propagation is simulated at times after steady-state 
conditions are reached but much earlier indicated by the experiments. 

The crack growth simulation results are shown in Fig. 5.19. The crack 
increment size for this study was approximately 0,01 mm which is nominally 2. A 
percent of the crack-tip element size. It can be seen that only one mesh shift 
(i.e. , crack growth step) was modeled. It is clear from this figure that the 
time it takes for the specimen to return to steady-state is significantly less 
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than the time to the next crack growth increment (indicated by dashed lines). 
Therefore^ the initial portion of the crack growth histories of Fig. 5.18 are 
clearly occurring under essentially steady-state conditions and thus C* as 
well as (T-) are valid crack-tip parameters. Since an Increase in C* re- 

X C X 

suits in a more rapid return to steady-state conditions, the above conclusion 
will remain valid for the initial constant velocity portions of the curves of 
Fig. 5.18. 

When crack growth occurs so slowly that the crack-tip is essentially 
at steady-state, the crack-tip field does not depend on the history of the 
specimen. Therefore, assuming steady-state conditions continue to exist, it is 
possible to skip to the final stages of crack growth without modeling the 
intermediate crack growth. If it is found that crack growth is still slow 
enough for steady-state conditions to exist, then it seems reasonable to expect 
that the bevalvor at Intermelate crack lengths is also of a steady-state type. 

The following describes the results of this procedure when applied to the 
two double-edge-crack specimens. 

Analysis of Final Stage of Crack Growth 

To analyze the final stage of crack growth, the crack length is in- 
creased to 2.75 mm and the process of applying the load elastically and creeping 
to steady-state is repeated. Table 5.3 summarizes the computational aspects of 
this process. The convergence of and C* to their steady-state values is 

shown in Fig. 5.20. Having reached steady-state, the cracks are grown at the 
rate suggested by the last portion of the crack histories (Fig. 5.18) as shown 
in Fig. 5.20. The significant Increase in the frquency of mesh shifting (compared 
to that in Fig. 5.19) due to the velocity Increase makes the details of the 
curve difficult to distinguish in this figure. However, the time step size is 
such that six or seven steps occur between each crack growth increment. Unlike 
the strip problem, the values of (Tj^)^ and C* are clearly increasing during 
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this qrack propagation process. 

It Is necessary to determine whether this Increase in the crack-tip 
parameters is due to the crack-tip no longer being at steady-state conditions 
or whether it is due to the increase in crack length* This is accomplished 
by continuiiig the calculation without further crack extension. If the value 
of the parameters do not change significantly with time, this means the increase 
was largely due to the crack length Increase and that crack growth is still 
occurring under essentially steady-state conditions* Examination of the final 
portions of the curves of Fig* 5.20 shows that this is the case* 

Based on this analysis; it appears that the conclusions reached as a 
result of the strip calculations are still valid. Since, (i) the strip analy- 
ses are much less expensive than this analysis of the double-edge-crack geometry, 
(ii) the steady-state Qj for the atrip is easily obtained analytically and (lii) 
the crack-tip parameters do not depend on crack length for the strip geometry, 
it seems that similar studies for other materials and/or other temperatures could 
most effectively be accomplished through the use of the strip geometry* The 
need for such studies follows from the vast simplification of fracture analysis 
and prediction which results if crack growth occurs under steady-state conditions* 
More will be said about this point in the conclusions* 
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SECTION VI 


CONCLUSIONS 


f JOBuary of Result o 

A finite element model has been derived which Is generally applicable 
to viscoplastic material models. This model uses an Initial strain approach 
which reduces computation time spent on forming and decomposing stiffness 
matrices and also circumvents the problem of element Incompressiblity constraints 
Through special features, including a correction term in the finite element 
equation, this model provides for Improved adherence to the postulated constitu- 
tive behavior (as compared to the standard initial strain approach) and 
allows time steps which approach in size those used in tangent stiffness 
methods. The accuracy and efficiency of this model with eight-node isopara- 
metric elements and the quarter-point crack-tip element approach have been 
verified through several calculations for a compact specimen geometry and a 
strip geometry. Also, a method of simulating crack growth througli. shifting 
of the quarter-point singularity elements and periodic remeshlng has been 
described and demonstrated. 

It has been shown that despite the fact that C* characterizes the 
crack-tip fields under steady-state creep conditions, it does not have an 
energy or energy rate interpretation. A related path-independent integral 
parameter however, does have the energy rate interpretation commonly 



attributed to C|, Since experimenteilsts use this energy interpretation to 
correlate creep crack growth rates, it appears that opposed to C|) 

la gaining acceptance as a useful creep crack growth rate criterion . Further- 
more, does not rely on the existence of steady-state creep conditions 

and thus might be expected to be a valid criterion even If creep crack growth 
should occur at rates which preclude the existence of steady-state creep 
conditions at the crack-tip. 

A creep crack growth simulation for 304 stainless steel has shown that 
for realistic load levels and corresponding crack speeds the crack-tip field 
Is essentially at a steady-statr; creep condition. This means that for this 
material, the propagating crack-tip field is largely unaffected by the 
history of crack growth or the history of loading. This feature can greatly 
reduce the analysis required for predlcltng creep crack growth behavior 
in a component as can be seen from the following suggested methodology. 

We assume that the crack propagation speed •j- is related to (T, ) 

• Qt X CSS 

^ E 

(i.c., - through the power law suggested by experimental data [23, 

24 ] . 




( 6 . 1 ) 


Nwxt we determine (e.g., by steady-state creep finite element analysis) 
('f. ) ^ as a function of crack length. Because of the assumed steady-state 
crack-tip behavior, this can be accomplished by considering several discrete 
crack lengths and then fitting a curve. No crack growth simulation procedures 
are necessiiry. Combining (6.1) with this result provides the following 
relationship between time and crack length: 




a(t) [(TJ ) 
L,jm.. 


-0 


da 4* t , 


( 6 . 2 ) 


w!:ere a^ is the iniriai crack length and t^ is the time when crack growth Initiates. 
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The* only unknown quantity In (6.2) is tha initiation time t , If the initiation 
tlrae foe creep crack growth is assumed to be negligible (as might be suggested 
from the results of [9,10]) then (6,2) immediately provides the predicted crock 
growth history. 

Vitek [11] does not consider t^ to be negligible based on severjil 
experiments (compact and double-edge-crack specimens) on two CrMoV steels. 

Using a dislocation model he further concludes that a measure of crack opening 
displacement (COD) correlates well with the initiation of crack gr«'*vth in these 
experiments. If the same conclusion is valid for 304 stainless steel, then 
one can presumably predict t^ based on a transient finite element analysis of 
the initial flawed configuration and a critical value of COD. If initiation 
occurs long after steady-state conditions are reached, it is then reasonable 
to estlr.^'te t^ using the rate of COD obtained from a steady-state finite 
element solution. The use of (6.2) and of the critical COD concept has not 
been Investigated in this study. 

All of the creep calculations have used the constitutive law which is 
obtained by generalizing the Norton constitutive law to three dimensions. 

Whereas this law is a good representation of steady-state creep behavior, it 
does not, in general, represent the primary stage of creep. Future work should 
include a study of other creep constitutive laws (such as that of Bodner and 
Fartom [32]). Also, the present model is derived on the assumption that dis- 
placements are small and strains infinitesimal. The strains in the vicinity of 
the crack-tip for the present calculations with 304 stainless steel material 
properties arc on the order of 5-lOJS and therefore suggest that a finite 
strain formulation may be more appropriate. A study should be undertaken to 
examine this aspect of the model. 

As noted previously, the creep crack growth prediction methodology 
expressed in (6.1) and (6.2) has not been tested. A study to assess the 
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utility o£ this rnethodoiogy should theroCore be undertaken. This study 
should consider crack growth Initiation as well as crock propagation and 
should Include a range of load levels and several specimen geometries. If 
the methodology is found to be successful for constant applied loads, then 
the study slvould be extended to consider more general load histories. 



APPENDIX A 


Existence of L im its for Contour Integral Definitions 
This appendix discusses the existence of the various llmito which have 
been taken In defining (AT) , (T) and C’*'. In considering these limits, we 

“ C ' C ■" ~ 

make use of the generally accepted result (see [6] for example) that the 
strain energy density quantities W and W as well as the quantity W* behave 
as 1/r In the vicinity of the crack-tip. This is assumed to be valid for 
nonsteady as well as steady-state creep and also for the elastic state 
existing at t ■ 0. 

Based on the known asymptotic behavior at the crack-tip (i.e., the 
HRR fields) the limits of contour integrals for equations (2.11,13,14,25,30, 
40,44) can be written in the following form provided one takes as being a 
circular contour centered at the crack-tip. 

(■J)f(£,e)ede - I f(e)d9 (A.l) 

A 

The nonsingular function f(e,0) becomes equal to f(9) when the limit is taken 
anc reflects the asymptotic nature of the HRR fields. 

In the following we limit the discussion to symmetric problems involving 
only mode 1 crack-tip deformation. Further, we assume that crack surface trac- 
tions and body forces are identically zero. With these conditions, we need only 

consider (AT.) and we can therefore rewrite (2.11) 

1 c 

^^'^1^ “ £-K) (A.2) 
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Now ('Ohiilvlot t*ho l imit of tho V^-V . IntogCiii of (A.2)i Inapcction of thia 


iivtogtiil ahowa that it can bo put into the form: 


1* **ki\ f “^r i|(i',o)tiU*do 

' ./v’-v 


U.3) 


t » 

R It 


I' J (^)6(t»0)iU'dv 


where V* i» » ymnll eircular region eentorod at the crack=tip and C is thn 


liUegrol over the region « V*. The function Is; a nonsingular 

ium ( f.otv tjhieh becomea g(h) In the limit aa r goes to Koro, whoro g(h) ia 
known In terma of the HRR fields# Upon a first inspection of (A, 3) one is 
t empted to eoiu'lude that* the limit does not exist since th® integrand has 
u non“integra!vle singularity at r *• 0. If, however ^ we look at the right 
ofiuallty of (A. 2) it is seen that this oonclusion results in a contradiction. 


Unee we have shown that the limit of the integral on f. does exist (and 


therefore (A. 2) requiros that the limit of the integral over V. - V, swat 
exlal), A ro~inspectlon of (A. 3) shows that the only way for this apparent 
eontraetlon to be resolved Is if tbn g(r,h) of (A. 3) has the following 
property. 


1! 11 

J ^(r,(V)dd 


(A.4X 


If function g(vl) la known explicitly for the linear elastic caao and thero' 
foro (A# 4) can be directly verified. For the HRR field g(d) is not known 

expllcltl.y and therefore (A. 4) can only be verified numerically. 
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For Inf inlt«8lmal strain, nonlinear elasticity, the following relation 
provldos an alteniatlvt,' to verifying (A. 4) Uireclty 

Du 

. „„ HS — / n Anr 

J 1,1 


/ 


Da, . 

— Lii. ^ ^ tjy 

Dx. i,j 


j ” ^“4 r 

”4 ^*^4 1 dS - I n 

I'.., J ’■'k J. 


23« 


Aa** dS (Aa5) 

J ij DXj^ 


t 

The relation (A# 5) (which assumes zero crack surface tractions and no body 
forces) illustrates that this volume integral of type (A. 3) can be expressed 
in terms of the contour integral of type (A.l). The relation (A. 5) can be 
verified througti the divergence theorem, the linear momentum balance condition 
and the following Identities: 


\ 
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Nuinei‘lcttl OiCference Metwoen (T, ) «t*d Cl? 

1 esa 1 

the pvu'poK© oC this appencUx la to give aome oxamplea to illuatrate 
the mimerteal cliff orencu hetwoen (T.) and C* «a given by (2,48). Ualng 

X *1 

(2,5(l)j (2,49a) and (2.51), we have 
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'fhe values tabulated in Table 8.1 weee computed approximately from values of 


1 and plots of (d) glvon in [6] and slcould bo viewed accordingly, 

eq 


Table 8. I Comparison of (T,) _ and C* 

Wane ^^raln Plane ^ Stress 
n * 3 n » 13 n » 3, n *■ 


(T, ) 


I csit 




0,98 1,00 


i.ll 1.14 


It Ah seen that for the range of n commonly encountered, (f.) and Ct are 

I CSiS X 


mimtjticany very slwAlar for plane strain but differ s,lgiiiftcantiy for plane 
stress. 


APPENDIX C 


Nui ngt’lcal Methods for Evaluation o| Contour Intagralti 
The numerical procediires for evaluating as defined by (2.27), 
as defined by (2.14) and CJ as defined by (2,31) are '•escribed in the following. 
General Procedures 

In studying the contour Integral paths indicated in the finite element 
meshes of Pigs. 5. 2, 3, 5 (dashed lines) it is sean that the paths always pass 
thruugli the centers of elements as opposed to along their edges. This procedure 
has been adopted so as to benefit from the presumably more accurate solution 
within tlie elements. Each element contour is divided into two segments with 
the integration being accomplished by two point Gaussian quadrature. All of 
the, Integrations are performed in the element local coordinates. 

The J I -Integral for linear Elastic Analyses 


The contour integral portion, of (2.27) Involves the stresses, » and 
the displacement derivatives, 8u^/9x^. Both of these quantities can be evalu- 
ated at tlie rctiuired Gauss points through the element nodal displacements and 
simple manipulations with element matrices. In the current study, is only 
considered as a parameter for linear elastic material behavior and therefore 

U » (l/2)a^jC^j. 

The C^- Integral 


The GJ Integral of (2.31) consists only of a contour integral. The 
of (2.31) is evaluated using (2.43). The gradient rates are approximated by 

OXj^ At 
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and tlu'ireforu nro average rataB for tl'o incremont as o|Vj)Qsed to the rates at 
tiu' end of the increment, ITie contour integration procedure for Cj In as 
described above and uses two point (luadrature for each element segment. Whereas 
stresses are easily computed at the riui/ilred contour Gauss points in the elastic 
case* tlve stresses must be computed incrementally in creep analyses and there- 
fore stress Information must be stored for each contour integration point 
unless nonstandard clement interpolations are used. Xn the present study, 
the stresses at the contour Gauss points are Interpolated from the 2x2 element 
Gauss points through bilinear Lagrangian interpolation (in local coordinates), 
thus eliminating the, need for additional storage. 


The ( Tj)^ -Integral 




In the evaluation of (2.14) It Is understood that are the stresses 
at the beginning of chc time .increment being considered. The procedures for 
evaluating the contour integral portion of (2.14) are the same as '«sed in evalu- 
ating G*. The Incrmental stress-work density, 4W, is computed from 

iW - (ty + I 

'file stress derivative appearing in the area integral of (2.14) is evaluated 
based on the 2x2 element Causa mint values and the as.sumptlon that the stresses 
arc dlstrib\itcd bilinearly with respect to element local coordinates. Elements 
which are entirely within are Integrated with the usual 2x2 Gauss quadrature. 
Elements which are only partially within V have each applicable quadrant In- 
tug rated by one point Gauss quadrature. 


APPENDIX D 


Slroultttlon of Crack Extenalon 

Modeling the propagation of a crack using the finite element method 
requires some special procedure for representing the creation of new crack 
surface. A common procedure is to relax the nodal forces at the crack-tip 
node, thus in effect allowing the crack to extend to the next n&de along its 
path of propagation. This relaxation process can be accomplished in one 
time step but usually is allowed to extend over several time steps due to the 
large change in nodal forces which is inherent in the process. The major attrac- 
tion of this node- release procedure is its simplicity. There are two draw- 
backs of this procedure which resulted in an alternate procedure being adopted 
in this study. The first is that the Increment in crack growth is directly 
determined by the nodal spacing in the mesh^ therefore restricting the flex- 
ibility one has in selecting a time step size, the mesh size and /or the 
number of nodal force relaxation steps. The second and perhaps more important 
drawback, is that the method is not adaptable to models which use crack-tip 
singularity elements. 

A typical mesh in the vicinity of a crack-tip is shown In Fig. D.l. The 
region A represents the region being modeled by singular crack-tip elements 
which In the present case remain centered on the crack-tip. The Type B 
elements are eight-noded isoparametric elements which distort so that the 
region A can remain centered un the crack-tip. The sequence of element configu- 
rations in Fig. D.l illustrates the shift ing/rcmeshlng procedure used in [53, 

54] and adopted here. The region A is moved by shifting nodes without altering 
element connectivity until the Type B element ahead of the region A becomes 
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overly dlstored. At this point, the elements In the vlelnity are redellned 
so that iurther shifting is possible. It can be seen that this procedure allows 
the increment in crack length to be arbitrarily small and does not involve 
release of nodes in the same sense as for the previously described node- 
release procedure. 

The added flexibility afforded by this shifting procedure does require 
some additional work. I’or example, In the creep crack growth application, nodal 
displacements and element integration point stresses are Interpolated. The 
method of interpolation which is employed in this procedure is discussed next. 

We consider chat the solution at time t, has been obtained and we now 
must find the solution at time t^. During the Interval 

grown by an amount ha. Since the crack growth simulation procedure requires 
that nodes be shifted, and since the solution at tj^ must be represented In terms 
of nodal and Causs point quantities for the shifted mesh, it is necessary to 
submit the affected nodes and Gauss point to an interpolation fitting pro- 
cedure , 

The simplest lnterpo.'i;,tlon procedure for nodal displacements and the one 
used In [53,54] as well as for calculations in the present study is one which 
directly uses the element shape functions. In this method, the nodal positions 
for the mesh at are located in the mesh at t^^. Knowing which element 
of the mesh at encompasses this new node position allows the ImmediaCe 
calculation of displacements by use of the element shape functions and the 
nodal quantities for the mesh at While this is a consistent procedure 

for Uansfurrlng the solution at to the mesh at t^, it should be understood 
that the transfer cannot be perfect. That this must be the case can be seen 
by considering that spatial derivatives of displacements, etc., are not continu- 
ous across element boundaries. Since the element boundaries change position 
during the shifting process, points which had continuous derivatives at t^^ will 
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l\avt’ diacontlnous derivatives In the meeh at and vice-versa. 

When the mesh la shifted, the elemeirit Gauss pulnts are also shifted*, this 
means the Causa points represent different material points before and after the 
shift « In order that the new Gauss point stresses accurately reprosent the cur- 
rent stress state, it is necessary to Interpolate stresses for the new Gauss 
point locations using the old Gauss point values and locations. The procedure 
for doing this is to assume the element stresses are distributed blllnearly with 
respect to the clement local coordinates. Then it becomes possible to use 
bilinear hagranglan interpolation polynomials and the 2x2 element Gauss point 
streases to interpolate within each clement. For all creep crack growth 
calcuXations in tivis study, the crack growth increment sizes were chosen small 
enough that the new Gauss point stresses for each shifted element were always 
the result of interpolation within that same eljment* 


APPENDIX a 


Analytical Evaluation of C* tor the Strip Problem 


This appondix briefly outlines the anlaytlcal evaluation of C| for tvQ 
infinite strip problems and then summarizes the results in tabular form. 

The first strip problem is that which is illustrated in Pig. 5.3. We shall 
refer to this problem as Case A. The second problem, or Case is similar 
fo Case A in e.ve ry respect except the top and bottom edges of the strip are 
"clamped'* rather than *'on rollers". Those boundary conditions are summarized 
as follows: 

Case A: Uy(x,h) ** -u^Cx.-h) «" t (E.l) 

tj^yCx.h) Tj^y(x,-h) » 0 

Case B: Uy(x.h) • *Uy(x,-h) • 5 (E.2) 

OjjCx.h) “ u^(x,-h) - 0 

The crack surfaces are traction-free in both cases. 

We can select a C|-lntegral contour which allows C^ to be evaluated 
quite easily.^ Consider a contour of rectangular chape which coincides with the 
top and bottom edges of the strip, extends far enough ahead of the crack-tip 
so as to be in a steady-state stress field which is unaffected by the presence 
of the crack-tip. and extends far enough behind the crack-tip so as to be In 
stress-free material. We now evaluate C*. as defined by (2.31). through the 
use of this contour. It con be seen that for both Case A and Case B. the 
horizontal portions of the contour at y •• h do not contribute to the 

\his procedure parallels chat used by Rice [SO] for the evaluation of J. in 
a similar elastic strip problem. 
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Ituegiul, nor dotni dm portion In the stress-tree material. At the vcrltcal 
portion of the contoer ahead of the crack«tip, the only non-zero term is that 
involving W^. Therefori*, It is seen that for both Case A and Case B we have 

Ct « 2WAh (E.3) 

I ^ 

where W* implies W>'' existing far ahead of the crack-tip. Using the boundary 
conditions (H.l) and (G.2) and the assumption of steady-state conditions, it 
is possible to evaluate the remote steady-state stresses, . Using (2.43) re 
suits In W* and* thus Gi*. The results of this exercise are summarized in Table 

iU X 

The corresponding linear elastic strip problem which is obtained by re- 
placing the displacement rate boundary conditions by the corresponding dln- 
placeraent boundary onditions has been treated in a similar manner. These re- 
sults are also given in Table E.l. 
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Table E 

.1 Analytical Solutions for 
Infinite Strip problem 

the 



Steady -State Creep 

linear Elasticity 

Case A 


iiKT 


1 o 

stress 

1 1 

1 

piano 

strain 

/2Vntl 2 

\s) S 

1 

1-v^ 

Case B 




plane 

stress 

I 2\n+t 2 

\B) 3 

1 

l-y2 

plane 

strain 

# # 

1-v 

(l+v)(l-2v ) 


This case does not have a steady-state solution since the 
boundary conditions require a volumetric strain rate. 
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K. Kathlresan), 280 pages. 

(15) "Fatigue Crack Growth In Modes 1 and XI Spectrum Loading" NASA- 
CR-78-123. (with M. Nakagakl) Oct. 1978, 98 pages. 

( 16 ) "Noz-Flaw: A Computer Program for Direct Evaluation of K-Pactors 

for Pressure-Vessel Nozzle Corner Cracks", (with R. Bar;i<>, J.W. 
Bryson, K. Kathlresan) NUREC/CR-18A3/ORNL/NUREG/CSD/TM-18, (Prep, 
for U.S. Nuclear Regulartory Comm, by Oak Ridge Natl. Labs), Nov. 
1980, 52 pages. 

(17) "Boundary Element Methods (BEM) and Combination of BEM-FEM" 

(with J.J. Grannell) GXT-CACM-SNA-79-16, 84 pages, 1979. 


Papers Presented at International and National Conferences 

(1) ''Nonlinear Free Oscillation of Shells", Regional Meeting of 
Society of Industrial and Applied Mathematics . November 1971, 

Ellensburg, WA 

(2) "Nonlinear Oscillations of a Hinged Beam Including Nonlinear Inertia 
Effects", ASME Joint National and Western Applied Mechanics Con- 
ference . June 26-29, 1972, University of California, San Diego. 

(3) "influence of Large Amplitudes and Boundary Conditions on the Super- 
sonic Flutter of a Cylindrical Shell", presented at the 4th 
Canadian Congress of Applied Mechanics . May 28-June 1, 1973, 

Montreal, Canada. 

(4) Mechanics of Brain Tissue Fragility", ASME Special Symposium on 
Biomechanicst National Summer Conference on Applied Mechanics . 

June "20-22, 1973, Atlanta, GA. 

(5) "Nonlinear Flutter of a Cylindrical Shell", 7th SECTAM . Catholic 
University of America, March 1974. 

(6) "Application of an Assumed Displacement Hybrid Finite Element 
Model to Two-Dimensional Problems in Fracture Mechanics", 

AlAA/ASME/SAE 15th SDM Specialist Conference . Las Vegas, NV, 

April 1974. 

(7) "An Assumed Displacement Hybrid Finite Element Model for Linear 
Fracture Mechanics", presented at the 7th U.S. National Congress 
of Applied Mechanics . Boulder, CO, June 1974. 

(8) "Finite Element Program for Fracture Mechanics Analysis of 
Composite Materials", presented at ASTM Symposium on Fracture of 
High Modulus Fibers and their Composites. Natioal Bureau of 
Standards . Gaithersburg, MD, Sept. 25, 1974 (invited). 

(9) "Stress Intensity Factors of Cra<?S<.ed Orthotropic Plates", 

Conference on the Fundamental Aspects of the Deformation an d 
Fracture of Composite Materials . Bat telle Seattle Research 
Center, Seattle, Feb. 22-24, 1975. 

(10) "Brain-Tissue Fragility: A Finite Strain Analysis by a Hybrid 

Finite Element Method", A merican Society of Mechanical 
Engineers: Applied Mechar.ics Western Conference . University of 

Hawaii, March 25-27, 1975. 

(11) "Rotationally Symmetric Bending of Orthotropic Conical Sheila: 

Transverse Shear and Couple Stress-Stress Couple Effects", 

14th Midwestern Mechanics Conference, University of Oklahoma, 

Norman OK, March 23-25. 

(12) "Three-Dimensional Cracked Elements", AFRPL Edwards Air Force Base, 

CA, Contract Research Review, Caltech, Pasadena, May 8-9, 1975 (invited) 


(13) "Boundary Integral Equation Eocmulatlon for Three-Dimensional 
Elasticity Ftoblems with Body Forces", 5th Canadian Congress of 
Applied Mechanics . New Brunswick, Canada, June 1975 (with A.L* Deak) 

(14) "Nonlinear Stress Analysis of Pneumatic Structures", 5 ^;h Canadian 
Congress of Applied Mechanics . New Brunswick, Canada, June 1975, 
(with A.L. Deak). 

(15) "Finite Element Approximation In Solid Mechanics", four lectures 
at the University of Tennessee Space Institute, Tullahoma, Nov. 
10-15, 1975. 

(16) "Large-Scale Yielding Fracture Mechanics", at the Committee E-24 
Meeting. 9th National Symposium on Fracture . University of 
Pittsburgh, Aug. 25-27, 1975 (invited). 

(17) "Analysis of Two-Dimensional Problems Involving Large-Scale 
Yielding", 12th Annual Meeting of the Society of Engineering 
Science , Austin, TX, Oct. 1975. 

( 18 ) "An Assumed Displacement Hybrid Finite Element Model for Three- 
Dimensional Linear Fracture Meclwnics Analysis", 12th Annual 
Society of Engineering Science Meeting . University of Texas at 
Austl^i, Oct, 1975* 

(19) "Hybrid Elements for 3-D Fracture", Specialist Work Shop on 3-D 
Fracture Analysis, organized by AFOSR/NASA/ERDA/DOT, Battelle 
Columbus Labs., Columbus, OH (invited). 

(20) "Post-Yield Analysis of a Three-Point -Bend Fracture Test Specimen", 
8th Southeastern Conference on Theoretical and Applied Mechanics . 

VPI and SU, Blacksburg, VA, April 20-30, 1976. 

(21) "J-Integral Estimation Procedures for Strain-Hardening Materials", 
at AIAA/ASME/SAE, Structures, Structural Dynamics, and Materials 
Specialist Conference . Valley Forge, PA, May 5-7, 1976. 

(22) "On a 3-D Singularity Element for Computation of Mixed Mode Stress 
Intensities", 13th Annual Society of Engineering Science Meeting . 
NASA-Langley, Hampton, VA, Nov. 1976. 

(23) "Fracture Analysis Under Large-Scale Plastic Yielding Conditions", 
10th U.S. National Conference in Fracture . American Society for 
Testing and Materials, Aug. 20-25, 1976, Philadelphia, PA. 

(24) "On Hybrid Stress Analysis of Laminated Shells by the Hybrid 
Stress Finite Element Model", presented at International Con- 
ference on Computational Methods in Nonlinear Mechanics , Univer- 
sity of Texas at Austin, Sept. 1974. 

(25) "On Hybrid Stress and Hybrid Displacement Models in Solid and 
Fracture Mechanics", AICA International Symposium on Computer 
Methods for Partial Differential Equations . Lehigh University, 
Bethlehem, PA, June 1975 (invited paper). 
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(2b) "Finite Element ^inalysis of Cracks Between Dissimilar Media", 

NATO Advanced Study Institute on Continuum Meclmnlcs Aspects of 
Geodynarolcs and Rock Fracture Mechanics, Reykjavik, Iceland, 

Aug. 11-20, 1974 (invited) . 

(27) "Finite Element -Perturbation Analysis on Nonlinear Dynamic Response 
of Elastic Contlnuua", invited presentation at the 1974 International 
Conference on Finite Element Methods in Engineering . Sydney, Aus- 
tralia, Sept. 1974, 

(28) "Three-Dimensional Linear Fracture Mechanics: Analysis by a 

Displacement Hybrid Finite Element Model", invited presentation 
of the 3rd International Conference on Structural Mechanics In 
Reactor Technology . University of London, Sept. 19 "^5, (with 

K, Kathlresan and A.S. Kobayashl) , 

(29) "Stress Analysis of Cracks In Elosto-Plastic Range", (with M. 

Nakagakl) , 4th Quadrennial International Conference on Fracture , 
University of Waterloo, Ontario, Canada, June 1977. 

(30) "Stress Intensity Factors for Surface Flews in Pressuclzed Cylinders", 
(with K. Kiithiresan) , 3rd International Congress on Pressure Vessel 
Technology, Tokyo, Japan, April 1977, 

(31) "Fracture Initiation in Plane Ductile Fracture Problems", 3rd In- 
ternational ConRress on Pressure Vessel Technology. Tokyo, Japan, 

April 1977, 

(32) "On Hybrid Techniques for Fracture Analysis", (with P. Tong), 
International Conference on Fracture Mechanics and Technology . 

Hong Kong, March 1977 (tnvitedy. 

(33) "On the Formulation and Application of Rational Numerical Methods 
for Problems with Nonremovable Singularities", International Sym - 
posium on Innovative Numerical Methods in Engineering Science . 

Paris, France, May 1977 (invited). 

(34) "Edge Function Method for Three-Dimensional Stress Analysis", 

(with P.M. Quinlan and J.E. Fitzgerald), at International Sym- 
posium on Innovative Numerical Methods In Engineering Science . 

Paris, France, May 1977 (Invited). 

(35) "Fracture Analysis of Structures Under Combined Mode Loading", 

2nd ASCE Engineering Mechanics Specialist Conference , Rayleigh, 

N.C., May 1977, (invited). 

(36) "Stress Analysis of Cracks in Elasto-Plastic Range", 4th Inter - 
national Conference on Fracture . University of Waterloo, Canada, 

June 1977. 

(37) "Outer Surface Flaws in Pressure Vessels, 4th International 
Conference on Structural Mechanics in Reactor T e chnology , San 
Francisco, CA, August 1977* 


(38) ''Through Flaws in Flatus in Bending", 4th International Con- 
ference on Structural Mechanics in Reactor Technology. San 
Francisco, CA, Aug. 1977. 

(39) "Hybrid Finite Element Models in Nonlinear Solid Mechanics", 
International Conference on Finite Elements in Nonlinear Solid 
and Structural Mechanics . Cello, Norway, August 30-Septembcr 1, 

1977 (Invited). 

(40) "Analysis of Stable Crack Growth in Ductile Materials", 9th 
SAMPE National Conference . Atlanta, GA, Oct. 1977. (Invited). 

(41) "Surface Flaws in Plates", I4th Annual Meeting of Society of 
Engineering Science . Lehigh University, Nov. 1977. 

(42) "A Finite Element Analysis of Stable Crack Growth - l", ASTM 
National Symposium on Elastic-Plastic Fracture . Atlanta, GA, 

Nov. 1977. 

(43) "Hybrid Finite Element Models in Linear and Nonlinear Fracture", 
International Conference on Numerical Methods in Fracture . 

Swansea, U.K., Jan., 1978, (invited). 

(44) "Bi-Material Fracture", 1978 Joint AFOSR/AFRPL Rocket Propulsion 
Research Meetin g, Edwards Air Force Base, CA, March, 1978 (Invited). 

(45) "Elastic-Plastic Analysis of 3-D Cracks", 1978 Joint AFOSR/AFRPL 
Rocket Propulsion Research Conference . Edwards Air Force Base, 

CA, March 1978, (Invited). 

(46) "Stress Analysis of Typical Flaws in Aerospace Structural Components" 
19th AIAA/ASME Structures, Structural Dynamics and MateriaU Con- 
ference . Bethesda, MD, March 1978. 

(47) "Edge-Function Method for 3-D Elasticity", 8th U.S. National 
Congress (Quadrennial) on Applied Mechanics , UCLA, Los Angeles, 

June 1978. 

(48) "Boundary-Discretization Method Using Edge Functions", International 
Conference on Recent Advances in Boundary Element Methods . Univer- 
sity Southampton, U.K., July 1978 (Invited). 

(49) "Complementary Energy Principles and Finite Strain Problems", 
Symposium of the International Union of Theoretical and Applied 
Mechanics on Variational Methods ; Evanston, IL, Sept. 78, (Invited). 

(50) "Numerical Modelling of Nonlinear Behavior of Soft Biological 
Materials", International Conference on Applied Numerical 
Modelling . University of Madrid, Spain, Sept. 78, (invited). 

(51) "On the Use of Stress Functions and Asymptotic Solutions in Solid 
and Structural Mechanics", Symposium on Future Trends in Computer - 
ized Analysis of Structures *s . Washington, DC, Nov. 78, (invited). 
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(52) "An Efficient Assumed Stress Finite EXement Method for Analysis 
ov Angle-Ply Laminates", 15th Meeting of Society of Engineer Ing 
Sc lence , Unlv. of Florida, DecT * 78 , 

(53) "Finite Elasticity Solutions Using Hybrid Finite Elements Based 
on a Complementary Energy Principle", ASME Winter Annual Meeting . 

San Francisco, Dec. *78, 

(54) "Influence of Flaw Shapes on Stress Intensity Factors for Beltline 
Cracks", National Congress on Pressure Vessel Technology . San 
Francisco, June 1979, 

(55) "Finite Elasticity Solutions Using Hybrid Finite Elements Based 

on a Complementary Energy Principle II. Incompressible Materials", 
1979 Joint Applied Mechanics . Fluids Conference, Niagara Falls, 

New York, June 1979. 

(56) "Finite Element Methods for Finite Strain Plasticity Problems In 
Metal forming". International Conference on Computational M e thods 
In Nonlinear Mechanics . (sponsored by NSF) Austin, TX, March 1979, 
(Invited) . 

(57) "Analytical Modelling of Surface Flaws", Fracture Research Sym- 
posium, SESA, Annual Meeting . San Francisco, June 1979, (invited). 

(58) "Computational Methods for Engineering Fracture Analyses", Inter- 
national Conference on Fracture Mechanics In Engineering. Ban ga lore , 
India, March 1979, (InvltedV. 

(59) "Nozzle-Vessel-Intersection Cracks Under Thermal Shock", Inter- 
national Conference on Structural Mechanics in Reactor Technology . 

W. Berlin, Aug. 1979, (Invited). 

(60) "Hybrid Finite Element Methods for 3-D and Nonlinear Fracture Prob- 
l®™s". En gineering Applications of the Finite Element Methods . 

Det Norsks Veritas, Havlk, Norway, June 1979, (invited). 

(61) "Complementary Energy and Finite Strain Plasticity", Advances in 
Theory and Practice of Finite Element Methods , Centennial Celebration 
of Chalmers University of Technology, Goteberg, Sweden, Aug. 1979, 
(invited) , 

(62) "Statlc/Dynamlc Analysis of Crack Propagation", 3rd ASCE Engineering 
Mechanics Specialty Conf . , Unlv. of Texas at Austin, Sept. 1979, 
(invited) . 

(63) "Finite Strain Plasticity Computations", Society of Engineering 
Science Annual Meeting . Northwestern Unlv., Evanston, IL, Sept. 

1979 (Invited). 

(64) "Selection of Finite Element Bases", Specialist Workshop on Finite 
Elements . Washington Unlv., St. Louis, MO, Nov. 1979, (Invited). 


(65) "Numerical Modeling of Nonlinear and Dynamic Crack Propagation", 
Symp. on Nonlinear and Dynamic Fracture . ASME WAM, Dec, 79, 
(invited). 

(66) "Stress Analysis of Holes in Composite Laminates", 21st AlAA 
Structures. Structural Dynamics and Materials Conference . Seattle, 
WA, May 1980. 

(67) "Dynamic Propagation of a Central Crack in a Finite Panel", 

Int. Conf. on Analytical and Experimental Fracture Mechanics . 

Rome, Italyj June 1980, (invited). 

(68) "Use of Stress Functions and Asymptotic Solutions in FE^f Analysis 
of Continuua", Symp. on New Concepts in FEM. 1981 Summer Annual 
Mechanics Meeting. Boulder. CO, June 1981, (invited). 

(69) "Recent Studies in Hybrid FEM for Solids and Fluids", Society 
of Engineering Scl. Meeting . Atlanta, GA, Dec. 1979, (invited). 

(70) "An Embedded Elliptical Flaw Subject to Arbitrary Loading, in an 
Infinite Medium", 15th Int. Congress on Theoretical and Applied 
Mechanics. lUTAM , Univ. of Toronto, Aug. 1979. 

(71) "Edge-Function Method for Burled Cracks", 2nd Int. Symp. on 
Innovative Numerical Analysis in Engineering . Montreal, Canada, 
June 1980, (invited). 

(72) An Analy>ais of and Some Observations on Dynamic Fracture in an 
Impact Specimen", 1981 Pressure Vessels and Piping Conf. . Denver, 
CO, June 1981, (invited). 

(73) "Finite Defr 'cmation Analysis of Shells, A Complementary Energy- 
Hybrid Approach", Symp. on Nonlinear Finite Element Analysis 

of Shells . ASME, WAM, Nov. 81, (Invited). 

(74) "Recent Studies on Dynamic, Inelastic, and 3-D Fracture Analyeia", 
U.S. -Japan Seminar on Damage Tolerance Evaluation, Honolulu, 
Hawaii, Dec. 81, (invited). 


Invito d Seminars and CQjIoqula In U.S, and Abroad 

(1) "PiniCti RlRtnent: Analysis of S|iolls"i Department of Aerospace 
Knglncering, University of Maryland, July 1969* 

( 2 ) "Analysis of Large Amplitude Elastoplastic Dynamics of Shalls", 

May 1970. 

(3) "Philosophical Implications of the Theory of Relativity", 

December 1970. 

( 4 ) "Noivlinear Oscillations in Certain Elastic Systems", Indian 
Institute of Technology, Kanpur, September 1971. 

( 5 ) "Recent Developments in Finite-Element Theory", Department 
of Aeronautical Engineering, Indian Institute of Science, 

BungaLorc, September 1971. 

(6) "Ho.Ilcopter Ground Resonance", Department of Aeronautics and 
Astronautics, University of Washington, October 1971, 

(7) "What an Applied Heclianlclan Can Do in Medicine", Center for 
Bioengineering, University of Washington, February 1972. 

(8) "Peristaltic Pumping", Center for Bioengineering, School of 
Medicine, University of Washington. 

(9) "Nonlinear Oscillation in a ClrcuJ.ar Cylindrical Shell", Depart- 
ment of Mechanical, Mechanics and Aerospace Engineering, Illinois 
Institute of Technology, May 1973. 

( 10 ) "Dynamic Stability of a Shell in Supersonic Flow", Department of 
Theoretical and Applied Mechanics, University of Illlnols-Urbana, 
Champaign, June 1973. 

(11) "Analysis of a Rolling Aircraft Tire", Research and Development 
Center, General Tire and Rubber Company, Akron, OH, Aug. 9, 1974. 

(12) "Computational Methods In Fracture Mechanics", National Aeronautical 
Laboratory, Bangalore, India, Aug. 21, 1974. 

(13) "Elastic-Plastic Fracture Mechanics", Department of Aeronautics, 
Indian Institute of Science, Bangalore, India, Aug. 21, 1974. 

( 14 ) "Perturbation Methods in Nonlinear Flutter", Department of Aero- 
nautics, Indian Institute of Science, Bangalore, India, Aug. 30, 
1974. 

( 15 ) "Singular Pecturbation Methods in Shell Theory", lecture at Short 
Course on Singular Perturbation in Methods, University of Tennessee 
Space Institute, Tullahoma, TN, Nov. 4-8, 1974. 


(16) "finite ELemtint Approximanlon In Solid Mechanics", four lectures 
at Short Course on Approximate Methods In Engineering and Applied 
Sciences, University of Tennessee Space Institute, Tullahoma, TN, 
Nov,, 10-14, 1975. 

(17) "Gomputatlonal Fracture Mechanics", College of Engineering, Boston 
University, Boston, MA, Feb. 24, 1976. 

(18) "Novel Methods for Analysis of Singularity Problems", University 
College of Cork, Ireland, Jan. 78, 

(19) "Approximate Methods of Analysis", 4 lectures. University of 
Tennessee Space Institute, Tullahoma, TN, March 78. 

(20) "Recent Developments in Finite Element Methods", 4 lectures, M.I.T. 
July 78. 

(21) "Numerical Methods in Fracture Mechanics", University College, 

Cork, Ireland, Aug. 1979. 

(22) "Finite Strain Inelasticity Analysis Via Complementary Energy", 
Center for Computational Mechanics, Washington University, 

St. Louis, MO, Jan 80. 

(23) "Selection of Finite Element Basis", Center for Computational 
Mechanics, Washington University, St. Louis, MO, Nov. 1980, 

(24) "Dynamic Fracture Analysis", Dept, Mechanical Engg,, University 
of Washington, Seattle, WA, May 1980. 

(25) "Path-Independent Integrals in Fracture Mechanics", National 
Tslng-Hua University, Hsinchu, Republic of China, Taiwan, 

May, 1981. 

(26) "Rate Complementary Energy Principles for Finite Element Analysis 
of Finite Elasticity", National Taiwan University, Taipei, 

Republic of China, Taiwan, May 1981. 

(27) "Path-Independent Integrals in Fracture Mechanics", Rutgers 
University, N.J., July 1981. 

(28) "Path-Independent Integrals in Fracture Mechanics", Naval 
Research Laboratory, Washington, DC, Oct. 1981. 


